!23456789a123456789b123456789c123456789d123456789e123456789f123456789g12 ! ! Copyright D.J. Jeffery, 2006jan01. ! ! rde.f for preparing radioactive decay data. ! ! References: ! \bibitem[Browne \& Firestone(1986)]{browne1986} ! Browne, E., \& Firestone, R. B. 1986, ! Table of Radioactive Isotopes (New York: ! John Wiley \& Sons, Inc.) ! % Radioacitve isotopes and decay chains and ! % radiations including all the complex stuff like ! % X-ray, atomic electron, internal bremstrahlung X-ray. ! ! \bibitem[Huo(1992)]{huo1992} Huo, J. 1992, ! Nuclear Data Sheets, 67, 523. ! % all about isobars of atomic mass 56, ! % especially Ni56 and Co56. ! % An update of Huo~\etal 1987. ! ! \bibitem[Jeffery(1998)]{jeffery1998} Jeffery, D. J. 1998, ! astro-ph/9811356 ! % the longer LS grey $\gamma$-ray transfer procedure ! % for supernovae. ! % Effectively rejected A\&SS, but just for being ! % prolix and (in the referee's opinion) worthless. ! % But the referee did not say anything was wrong. ! % So unpublished rather than preprint is best. ! % The later implies eventual publication which ! % may not be forthcoming. ! % There is a good review of radioactive decay energy ! % quantities in Section 5. ! ! ! LBNL Isotopes ! Project - LUNDS Universitet ! WWW Table of Radioactive Isotopes ! http://ie.lbl.gov/toi/nuclide.asp?iZA=270056 ! Co-56 ! ! \bibitem[Metcalf et al.(2004)]{metcalf2004} Metcalf, M., Reid, J., ! \&~Cohen, M. 2004, ! Fortran 95/2003 Explained ! (Oxford: Oxford University Press) (MRC) ! % Pretty good reference book. ! !----------------------------------------------------------------------- ! program rde use numerical use physical use astronomical implicit none ! integer, parameter :: nelevel=200 integer, parameter :: nxray=200 ! real (kind=nprecision) :: con_mev_erg real (kind=nprecision) :: con_thalf_tefold integer :: i,j,k,l,m,n integer :: idata integer :: ielevel real (kind=nprecision) :: small real (kind=nprecision) :: xmass_np_diff ! character :: aanuclide*10 ! Nuclide name. real (kind=nprecision) :: b_coef ! Energy generation B coefficient ! ! (erg/s/g). real (kind=nprecision) :: b_coef_lepton ! For betas and atomic el. real (kind=nprecision) :: b_coef_photon ! For gammas and X-rays. real (kind=nprecision) :: b_coef_sun ! Energy generation B coefficient ! ! (erg/s/mass_sun). real (kind=nprecision) :: betamax ! Maximum beta particle energy (MeV). real (kind=nprecision) :: betafrac ! Fraction of time there is a ! ! beta particle. real (kind=nprecision) :: beta_ke ! Average beta particle energy ! ! per beta (MeV). real (kind=nprecision) :: c_coef ! Energy generation C ! ! coefficient (erg/s/g). real (kind=nprecision) :: c_coef_lepton ! For betas and atomic el. real (kind=nprecision) :: c_coef_photon ! For gammas and X-rays. real (kind=nprecision) :: c_coef_sun ! Energy generation C ! ! coefficient (erg/s/mass_sun). real (kind=nprecision) :: elevel(10,nelevel) ! Energy level ! ! and beta/EC data. real (kind=nprecision) :: fatomic ! Fraction of energy in atomic ! ! electron KE. real (kind=nprecision) :: fbeta ! Fraction of energy in positron KE. real (kind=nprecision) :: fgamma ! Fraction of energy in gamma-rays. real (kind=nprecision) :: flepton ! Fraction of energy in leptons. real (kind=nprecision) :: fphoton ! Fraction of energy in photons. real (kind=nprecision) :: fxray ! Fraction of energy in X-rays. integer :: ibeta ! The sign of the beta particle. integer :: iz ! Atomic number Z. integer :: ia ! Atomic mass number A. integer :: ia_parent ! Atomic mass number A of parent. integer :: ixray ! Number of X-ray emissions. real (kind=nprecision) :: qval ! Total disintegration energy (MeV). real (kind=nprecision) :: qval_nu_omit ! The energy that doesn't escape ! ! as neutrinos (MeV). real (kind=nprecision) :: qvalerr ! Error in qval. real (kind=nprecision) :: qval_atomic ! The energy in atomic electrons. real (kind=nprecision) :: qval_beta ! The energy in betas. real (kind=nprecision) :: qval_gamma ! The energy in gammas. real (kind=nprecision) :: qval_gamma_direct ! The energy in gamma lines. real (kind=nprecision) :: qval_xray ! The energy in gammas. real (kind=nprecision) :: tefold ! e-folding time (days). real (kind=nprecision) :: tefold_parent ! e-folding time (days) ! ! of parent. real (kind=nprecision) :: thalf ! Half-life (days). real (kind=nprecision) :: thalferr ! Half-life error (days). real (kind=nprecision) :: thalf_parent ! Half-life (days) of parent. real (kind=nprecision) :: xray(2,nxray) ! X-ray data. real (kind=nprecision) :: zzfiction ! Just a place-holder variable. ! !----------------------------------------------------------------------- ! namelist/param/aanuclide & & ,b_coef & & ,b_coef_lepton & & ,b_coef_photon & & ,b_coef_sun & & ,betamax & & ,betafrac & & ,beta_ke & & ,c_coef & & ,c_coef_lepton & & ,c_coef_photon & & ,c_coef_sun & & ,fatomic & & ,fbeta & & ,fgamma & & ,flepton & & ,fphoton & & ,fxray & & ,ibeta & & ,ia & & ,ia_parent & & ,iz & & ,iz_parent & & ,qval & & ,qval_nu_omit & & ,qvalerr & & ,qval_atomic & & ,qval_beta & & ,qval_gamma & & ,qval_xray & & ,tefold & & ,tefold_parent & & ,thalf & & ,thalferr & & ,thalf_parent & & ,zzfiction ! namelist/param2/elevel,xray ! !----------------------------------------------------------------------- write(*,*) 'Before preliminary results calculated.' !----------------------------------------------------------------------- ! con_mev_erg=1.d+6*echarge*1.d+7 ! 1.d+6 for MeV to eV, 1.d+7 to ergs. con_thalf_tefold=1.d0/log(2.d0) ! Convert half-life to e-folding time. small=1.d-10 xmass_np_diff=x_neutron_mass_mev-proton_mass_mev ! write(*,*) write(*,*) 'con_mev_erg,con_thalf_tefold' write(*,*) con_mev_erg,con_thalf_tefold ! !----------------------------------------------------------------------- write(*,*) 'Before data-read and evaluation of the quantities.' !----------------------------------------------------------------------- ! open(unit=1,file='rde.dat',status='old',action='read') ! call readcopy(1,0,0,'',idata) i=0 do410 : do elevel(:,:)=0.d0 read(1,param,end=210) read(1,param2,end=210) i=i+1 write(*,*) write(*,*) 'Nuclide ',i ! tefold=thalf*con_thalf_tefold ! !----------------------------------------------------------------------- ! write(*,*) 'Before counting the energy levels and xray ', ! write(*,*) 'emissions (IB plus lines).' !----------------------------------------------------------------------- ! j=1 do420 : do j=1,nelevel if(elevel(1,j) .lt. small) exit do420 end do do420 ielevel=min(j-1,nelevel)! The min function is redundant: j-1 would do. do430 : do j=1,nxray if(xray(1,j) .lt. small) exit do430 end do do430 ixray=min(j-1,nxray) ! The min function is redundant: j-1 would do. ! elevel(1,:ielevel)=elevel(1,:ielevel)*1.d-3 ! Conversion to MeV. elevel(3,:ielevel)=elevel(3,:ielevel)*1.d-2 ! Conversion to fraction. elevel(5,:ielevel)=elevel(5,:ielevel)*1.d-2 ! Conversion to fraction. qval_gamma_direct=sum( elevel(1,:ielevel) & & *(elevel(3,:ielevel)+elevel(5,:ielevel)) ) betafrac=sum(elevel(3,:ielevel)) qval_xray=sum( .001d0*xray(1,:ixray) & & *xray(2,:ixray)*.01d0) ! With conversion to MeV and ! ! fraction. ! !----------------------------------------------------------------------- ! ! The tricky part is the neutrino energy. We are assuming it escapes, ! and so we don't want to count it. ! But it is always included in the standard total disintegration energy. ! ! So we get the direct gamma fashion as above by counting up the energy ! That must come from decays from the excited child. ! ! We add to this, the mean beta kinetic energy per decay (not per ! beta particle) that we get from a table. ! ! We add to this, the pair annihilation energy if a beta+ is created. ! We assume that the positron is always annihilated. ! ! We neglect all atomic binding energies. See Enge-304--305. ! They are not always negligible, but they are at another order of ! accuracy than here. ! ! But we could account for them using the X-ray and Auger electron ! emission. But the Auger electron emission can be clumped with the ! internal conversion electron emission which competes with gamma-ray ! emission. Thus, one is in danger of double counting ! some energy both as gamma-rays and as internval conversion ! electrons. It's a bloody mess trying to sort it all out. ! There is also internal bremstrahlung X-ray emission. ! Better leave all this to anther day. ! But can get the data to do it in browne1986. It may be oline too. ! One should probably have xray and atelectron and photon and lepton ! categories. ! ! We also neglect nucleus recoil energy. It's probably pretty darn ! negligible since no one ever talks about it. ! ! Note in beta+ and EC, one goes from nucleus 1 and Z electrons in the ! universe to nucleus 2 and Z-1 electrons in the universe after annihilation. ! ! The total disintegration energy is atom 1 to atom 2 which also is ! Z to Z-1 electrons in the universe. So aside from binding eneryg ! everything is accurate. ! ! Now in beta- decay, one goes from nucleus 1 and Z electrons in the universe ! to nucleus 2 and Z+1 electrons in the universe. ! ! The total disintegration energy is atom 1 to atom 2 which also is ! Z to Z+1 electrons in the universe. So aside from binding eneryg ! everything is accurate. ! ! Recall the reactions: ! EC e**(-) + p**(+) --> n + neutrino ! beta+ decay p**(+) --> n + e**(+)+neutrino (Enge-302) ! beta- decay n --> p**(+) + e**(-)+antineutrino (Enge-303) !----------------------------------------------------------------------- ! qval_nu_omit=qval_gamma_direct & & +qval_beta & & +2.d0*emass_mev*betafrac*real(max(ibeta,0),nprecision) & & +qval_atomic & & +qval_xray fatomic=qval_atomic/qval_nu_omit ! Fraction of KE energy in atomic electrons. fbeta=qval_beta/qval_nu_omit ! Fraction of KE energy in positrons. flepton=fatomic+fbeta ! Fraction of KE energy in leptons fxray=qval_xray/qval_nu_omit ! Fraction of KE energy in x_rays. fgamma=1.d0-flepton-fxray ! Fraction of energy in gamma-rays fphoton=1.d0-flepton ! Fraction of energy in photons. qval_gamma=fgamma*qval_nu_omit if(betafrac .ge. small) then beta_ke=qval_beta/betafrac ! This is the mean beta KE per beta. else beta_ke=0.d0 end if ! c_coef=(qval_nu_omit*con_mev_erg) & & /(amu*real(ia,nprecision)*tefold*daysec) c_coef_lepton=c_coef*(fatomic+fbeta) c_coef_photon=c_coef*(fgamma+fxray) c_coef_sun=c_coef*xmass_sun if(ia_parent .ne. 0) then tefold_parent=thalf_parent*con_thalf_tefold b_coef=(qval_nu_omit*con_mev_erg) & & /(amu*real(ia,nprecision)*(tefold-tefold_parent)*daysec) else tefold_parent=0.d0 b_coef=0.d0 end if b_coef_lepton=b_coef*flepton b_coef_photon=b_coef*fphoton b_coef_sun=b_coef*xmass_sun ! write(*,param) ! end do do410 210 continue close(unit=1) ! end program rde ! !23456789a123456789b123456789c123456789d123456789e123456789f123456789g12 ! function decay(t) use numerical use physical use astronomical implicit none ! real (kind=nprecision) :: decay real (kind=nprecision), parameter :: ccco=6.70d9*xmass_sun ! Convert to per solar masses. real (kind=nprecision), parameter :: ccni=3.94d10*xmass_sun ! Convert to per solar masses. real (kind=nprecision), parameter :: fpe=3.20d-2 real (kind=nprecision), parameter :: fph=.968d0 ! real (kind=nprecision), parameter :: ftot=fpe+fph ! In principle 1 exactly. real (kind=nprecision), parameter :: ftot=1.d0 real (kind=nprecision), parameter :: gg=0.184641 real (kind=nprecision), parameter :: teco=111.48 real (kind=nprecision), parameter :: teni=8.767 ! real (kind=nprecision) :: decay real (kind=nprecision) :: t ! print*,'ccco,ccni' print*,ccco,ccni ! ! decay=6.31d43*exp(-t/8.8d0) + 1.43d43*exp(-t/111.d0) ! This is wrong in H2006. ! decay=ccni*( exp(-t/teni) + gg*(exp(-t/teco)-exp(-t/teni))*ftot ) ! end function decay