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      subroutine kap_get( &
 2         handle, species, chem_id, net_iso, xa, &
 3         logRho, logT, &
 4         lnfree_e, d_lnfree_e_dlnRho, d_lnfree_e_dlnT, &
 5         eta, d_eta_dlnRho, d_eta_dlnT , &
 6         kap_fracs, kap, dlnkap_dlnRho, dlnkap_dlnT, dlnkap_dxa, ierr)
 7         use chem_def, only: chem_isos
 8         use chem_lib, only: basic_composition_info
 9         use kap_def, only : kap_is_initialized, Kap_General_Info, num_kap_fracs, i_frac_Type2
10         use kap_eval, only : Get_kap_Results
11         ! INPUT
12         integer, intent(in) :: handle ! from alloc_kap_handle; in star, pass s% kap_handle
13         integer, intent(in) :: species
14         integer, pointer :: chem_id(:) ! maps species to chem id
15         integer, pointer :: net_iso(:) ! maps chem id to species number
16         real(dp), intent(in) :: xa(:) ! mass fractions
17         real(dp), intent(in) :: logRho ! density
18         real(dp), intent(in) :: logT ! temperature
19         real(dp), intent(in) :: lnfree_e, d_lnfree_e_dlnRho, d_lnfree_e_dlnT
20         ! free_e := total combined number per nucleon of free electrons and positrons
21         real(dp), intent(in)  :: eta, d_eta_dlnRho, d_eta_dlnT
22         ! eta := electron degeneracy parameter from eos
23
24         ! OUTPUT
25         real(dp), intent(out) :: kap_fracs(num_kap_fracs)
26         real(dp), intent(out) :: kap ! opacity
27         real(dp), intent(out) :: dlnkap_dlnRho ! partial derivative at constant T
28         real(dp), intent(out) :: dlnkap_dlnT   ! partial derivative at constant Rho
29         real(dp), intent(out) :: dlnkap_dxa(:) ! partial derivative w.r.t. species
30         integer, intent(out) :: ierr ! 0 means AOK.
31
32         type (Kap_General_Info), pointer :: rq
33
34         real(dp) :: X, Y, Z, abar, zbar, z2bar, z53bar, ye, mass_correction, sumx
35         real(dp) :: XC, XN, XO, XNe
36         integer :: i, iz
37
38         ierr = 0
39         if (.not. kap_is_initialized) then
40            ierr=-1
41            return
42         endif
43
44         call kap_ptr(handle,rq,ierr)
45         if (ierr /= 0) return
46
47         call basic_composition_info( &
48            species, chem_id, xa, X, Y, Z, &
49            abar, zbar, z2bar, z53bar, ye, mass_correction, sumx)
50
51         xc = 0; xn = 0; xo = 0; xne = 0
52         do i=1, species
53            iz = chem_isos% Z(chem_id(i))
54            select case(iz)
55            case (6)
56               xc = xc + xa(i)
57            case (7)
58               xn = xn + xa(i)
59            case (8)
60               xo = xo + xa(i)
61            case (10)
62               xne = xne + xa(i)
63            end select
64         end do
65
66         call Get_kap_Results( &
67            rq, zbar, X, Z, XC, XN, XO, XNe, logRho, logT, &
68            lnfree_e, d_lnfree_e_dlnRho, d_lnfree_e_dlnT, &
69            eta, d_eta_dlnRho, d_eta_dlnT, &
70            kap_fracs, kap, dlnkap_dlnRho, dlnkap_dlnT, ierr)
71
72         ! composition derivatives not implemented
73         dlnkap_dxa = 0
74
75      end subroutine kap_get