kap module InterfaceΒΆ

The primary entry point to the kap module is through the routine kap_get. Broadly, one provides the density, temperature and full composition information, as well as information about the free electron fraction and electron chemical potential (required for calculation of the Compton scattering opacity). Then the module returns the opacity and its temperature and density derivatives.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
      subroutine kap_get( &
         handle, species, chem_id, net_iso, xa, &
         logRho, logT, &
         lnfree_e, d_lnfree_e_dlnRho, d_lnfree_e_dlnT, &
         eta, d_eta_dlnRho, d_eta_dlnT , &
         kap_fracs, kap, dlnkap_dlnRho, dlnkap_dlnT, dlnkap_dxa, ierr)
         use chem_def, only: chem_isos
         use chem_lib, only: basic_composition_info
         use kap_def, only : kap_is_initialized, Kap_General_Info, num_kap_fracs, i_frac_Type2
         use kap_eval, only : Get_kap_Results
         ! INPUT
         integer, intent(in) :: handle ! from alloc_kap_handle; in star, pass s% kap_handle
         integer, intent(in) :: species
         integer, pointer :: chem_id(:) ! maps species to chem id
         integer, pointer :: net_iso(:) ! maps chem id to species number
         real(dp), intent(in) :: xa(:) ! mass fractions
         real(dp), intent(in) :: logRho ! density
         real(dp), intent(in) :: logT ! temperature
         real(dp), intent(in) :: lnfree_e, d_lnfree_e_dlnRho, d_lnfree_e_dlnT
         ! free_e := total combined number per nucleon of free electrons and positrons
         real(dp), intent(in)  :: eta, d_eta_dlnRho, d_eta_dlnT
         ! eta := electron degeneracy parameter from eos

         ! OUTPUT
         real(dp), intent(out) :: kap_fracs(num_kap_fracs)
         real(dp), intent(out) :: kap ! opacity
         real(dp), intent(out) :: dlnkap_dlnRho ! partial derivative at constant T
         real(dp), intent(out) :: dlnkap_dlnT   ! partial derivative at constant Rho
         real(dp), intent(out) :: dlnkap_dxa(:) ! partial derivative w.r.t. species
         integer, intent(out) :: ierr ! 0 means AOK.

         type (Kap_General_Info), pointer :: rq

         real(dp) :: X, Y, Z, abar, zbar, z2bar, z53bar, ye, mass_correction, sumx
         real(dp) :: XC, XN, XO, XNe
         integer :: i, iz

         ierr = 0
         if (.not. kap_is_initialized) then
            ierr=-1
            return
         endif

         call kap_ptr(handle,rq,ierr)
         if (ierr /= 0) return

         call basic_composition_info( &
            species, chem_id, xa, X, Y, Z, &
            abar, zbar, z2bar, z53bar, ye, mass_correction, sumx)

         xc = 0; xn = 0; xo = 0; xne = 0
         do i=1, species
            iz = chem_isos% Z(chem_id(i))
            select case(iz)
            case (6)
               xc = xc + xa(i)
            case (7)
               xn = xn + xa(i)
            case (8)
               xo = xo + xa(i)
            case (10)
               xne = xne + xa(i)
            end select
         end do

         call Get_kap_Results( &
            rq, zbar, X, Z, XC, XN, XO, XNe, logRho, logT, &
            lnfree_e, d_lnfree_e_dlnRho, d_lnfree_e_dlnT, &
            eta, d_eta_dlnRho, d_eta_dlnT, &
            kap_fracs, kap, dlnkap_dlnRho, dlnkap_dlnT, ierr)

         ! composition derivatives not implemented
         dlnkap_dxa = 0

      end subroutine kap_get