Line data    Source code

       1              : ! ***********************************************************************
       2              : !
       3              : !   Copyright (C) 2013  The MESA Team
       4              : !
       5              : !   This program is free software: you can redistribute it and/or modify
       6              : !   it under the terms of the GNU Lesser General Public License
       7              : !   as published by the Free Software Foundation,
       8              : !   either version 3 of the License, or (at your option) any later version.
       9              : !
      10              : !   This program is distributed in the hope that it will be useful,
      11              : !   but WITHOUT ANY WARRANTY; without even the implied warranty of
      12              : !   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
      13              : !   See the GNU Lesser General Public License for more details.
      14              : !
      15              : !   You should have received a copy of the GNU Lesser General Public License
      16              : !   along with this program. If not, see <https://www.gnu.org/licenses/>.
      17              : !
      18              : ! ***********************************************************************
      19              : 
      20              :       module neu
      21              : 
      22              :       use star_private_def
      23              :       use const_def, only: dp
      24              :       use utils_lib
      25              : 
      26              :       implicit none
      27              : 
      28              :       private
      29              :       public :: do_neu_for_cell, do_clear_neu_for_cell
      30              : 
      31              :       contains
      32              : 
      33        53310 :       subroutine do_clear_neu_for_cell(s,k,ierr)
      34              :          type (star_info), pointer :: s
      35              :          integer, intent(in) :: k
      36              :          integer, intent(out) :: ierr
      37        53310 :          ierr = 0
      38        53310 :          s% non_nuc_neu(k) = 0
      39        53310 :          s% d_nonnucneu_dlnd(k) = 0
      40        53310 :          s% d_nonnucneu_dlnT(k) = 0
      41        53310 :          s% nonnucneu_plas(k) = 0
      42        53310 :          s% nonnucneu_brem(k) = 0
      43        53310 :          s% nonnucneu_phot(k) = 0
      44        53310 :          s% nonnucneu_pair(k) = 0
      45        53310 :          s% nonnucneu_reco(k) = 0
      46        53310 :       end subroutine do_clear_neu_for_cell
      47              : 
      48              : 
      49        26684 :       subroutine do_neu_for_cell(s,k,ierr)
      50              :          use neu_def
      51              :          use neu_lib
      52              :          use const_def, only:ln10
      53              :          type (star_info), pointer :: s
      54              :          integer, intent(in) :: k
      55              :          integer, intent(out) :: ierr
      56              : 
      57       160104 :          real(dp) :: loss(num_neu_rvs)  ! total from all sources
      58       827204 :          real(dp) :: sources(num_neu_types, num_neu_rvs)
      59        26684 :          real(dp) :: log10_rho, log10_T
      60              :          real(dp), parameter :: log10_Tlim = 7.5d0
      61              :          logical :: flags(num_neu_types)  ! true if should include the type of loss
      62              : 
      63              :          include 'formats'
      64              : 
      65              :          ierr = 0
      66              : 
      67        26684 :          if (s% non_nuc_neu_factor <= 0d0) then
      68            0 :             call do_clear_neu_for_cell(s,k,ierr)
      69            0 :             return
      70              :          end if
      71              : 
      72       160104 :          flags = .true.
      73              :          !flags(reco_neu_type) = .false.
      74              : 
      75        26684 :          log10_rho = s% lnd(k)/ln10
      76        26684 :          log10_T = s% lnT(k)/ln10
      77              : 
      78        26684 :          if (s% use_other_neu) then
      79              :             call s% other_neu( &
      80              :                s% id, k, s% T(k), log10_T, s% rho(k), log10_rho, &
      81              :                s% abar(k), s% zbar(k), &
      82            0 :                log10_Tlim, flags, loss, sources, ierr)
      83              :          else
      84              :             call neu_get( &
      85              :                s% T(k), log10_T, s% rho(k), log10_rho, &
      86              :                s% abar(k), s% zbar(k), &
      87        26684 :                log10_Tlim, flags, loss, sources, ierr)
      88              :          end if
      89              : 
      90        26684 :          if (ierr /= 0) then
      91            0 :             if (s% report_ierr) then
      92            0 :                write(*,3) 'do_neu_for_cell: neu_get ierr', ierr, k
      93            0 :                write(*,1) 'T=', s% T(k)
      94            0 :                write(*,1) 'log10_T=', log10_T
      95            0 :                write(*,1) 'rho=', s% rho(k)
      96            0 :                write(*,1) 'log10_rho=', log10_rho
      97            0 :                write(*,1) 'abar', s% abar(k)
      98            0 :                write(*,1) 'zbar', s% zbar(k)
      99            0 :                write(*,'(A)')
     100            0 :                return
     101              :                stop
     102              :             end if
     103              :             return
     104              :          end if
     105              : 
     106        26684 :          if (s% non_nuc_neu_factor /= 1d0) loss(:) = loss(:)*s% non_nuc_neu_factor
     107        26684 :          s% non_nuc_neu(k) = loss(ineu)
     108        26684 :          s% d_nonnucneu_dlnd(k) = loss(idneu_dRho)*s% rho(k)
     109        26684 :          s% d_nonnucneu_dlnT(k) = loss(idneu_dT)*s% T(k)
     110              : 
     111        26684 :          s% nonnucneu_plas(k) = sources(plas_neu_type,ineu)
     112        26684 :          s% nonnucneu_brem(k) = sources(brem_neu_type,ineu)
     113        26684 :          s% nonnucneu_phot(k) = sources(phot_neu_type,ineu)
     114        26684 :          s% nonnucneu_pair(k) = sources(pair_neu_type,ineu)
     115        26684 :          s% nonnucneu_reco(k) = sources(reco_neu_type,ineu)
     116              : 
     117        26684 :          if (is_bad(s% non_nuc_neu(k))) then
     118            0 :             ierr = -1
     119            0 :             if (s% report_ierr) write(*,*) 'do_neu_for_cell ierr for cell', k
     120            0 :             if (s% stop_for_bad_nums) then
     121            0 :                write(*,2) 'bad s% non_nuc_neu(k)', k, s% non_nuc_neu(k)
     122            0 :                call mesa_error(__FILE__,__LINE__,'do_neu_for_cell')
     123              :             end if
     124            0 :             return
     125              :          end if
     126              : 
     127              :       end subroutine do_neu_for_cell
     128              : 
     129              :       end module neu