Hot Fixes to C90.02
These are minor changes that will be incorporated in the next revision of the code. These only apply to the current version of the code stored on this site. In the following, the corrections that should be made to version 90.02 are indicated by line numbers within specific subroutines.
Incorrect index on evaluation of derivative of photoionization heating. This did not affect results but did hinder convergence of the code, and could cause problems. Change line 106 of routine sumheat to the following: (change i-1 to i).
doi=1,limelm *hotfix 1997 Mar 18, sum should be to i, caught by Valentina Luridiana doj=1,i
level3 is not well protected against underflow A more bullet-proof version of that routine can be downloaded from here. This does not affect the results at all. so if only needed if the code crashes due to division by zero within level3.
Wrong stage of ionization multiplies HeI singlet cooling. Two changes should be made to line 388 of routine coolr. The heiii should be changed to heii and the multiplication by 4 should be removed. C90.03 will have this treatment improved.
2 (he1n(l)*he1bol(l,iu)*he1stat(iu)/he1stat(l)-he1n(iu))* 3 eden*heii* (1.384e-16* HdeltaT(iu,l) ) *hotfix 97 mar 4, above had been heiii, caught by Anthony Leonard *factor of 4 was also for He+, should be hydrogenic
Wrong temperature appears in calculation of ionic recombination on grain surfaces. The grain temperature appears, but the electron temperature should be there. Edit the file hmole.f and change line 365 to
* gionrc = 5.8e-13 / sqrt( tdlow) * drate * hden * *hot fix 97 feb 24, sqrt(tdlow) should be sqrte, caught by Jon Slavin gionrc = 5.8e-13 / sqrte * drate * hden * 1 grecon
Factor of electron density was missing from diffuse recombination rates of the heavy elements. This caused the effects of the fields to be underestimated for electron densities greater than one, and overestimated for electron densities less than one. In some very low density cases energy was not conserved, a problem caught by the code's internal sanity checks. Edit routine MakeRecomb and change lines 146 from
RecomRate(i) = eden*( radrec + rcota(i) ) * this is supposed to be the ground term of the recombination DiffHeavy(nelem,i ) = radrec / 10. xLyaHeavy(nelem,i) = radrec / 10.
RecomRate(i) = eden*( radrec + rcota(i) ) * this is supposed to be the ground term of the recombination DiffHeavy(nelem,i ) = eden*radrec / 10. xLyaHeavy(nelem,i) = eden*radrec / 10.
that is - add the factor of eden on two lines. 1997 Feb 17.
Escape probabilities were not computed on second and later iterations if the first iteration was fully ionized. This In routine MakeStatRT change line 161 from
* possible for line to mase, reason for sec test if( t(ipLnTauTot)*0.9-t(ipLnTauIn) .lt. 0. .and. 1 t(ipLnTauIn).gt.0. ) then
*change 97 jan 8, on iterations where the first iteration had NO species, *but it does exist on later iterations, added test for existence of ion if( t(ipLnTauTot)*0.9-t(ipLnTauIn) .lt. 0. .and. 1 t(ipLnTauIn).gt.0. .and. t(ipLnTauTot).ne.taumin ) then
It is also necessary to include the common block /taumin/; go to line 134 and add the two lines
common/DstWght/ DstWght , DstWghtSave real taumin , tlamin common/taumin/ taumin , tlamin
(just add the last two lines).
H charge transfer with highly ionized species Versions 90.02 and before did not include charge transfer recombination of highly ionized species with atomic hydrogen, because rates had not been computed. Alex Dalgarno has computed these rates (Ferland, Korista, Verner, and Dalgarno, ApJ Letters June 1 1997) and finds a large rate coefficient. This process can be easily added, and does effect results for photoionization by a hard continuum. To add charge transfer recombination for ionized species, edit the routine hctrecom, going to line 29. Change
* * return zero for highly ionized species if( ipIon.gt.4 ) then HCTRecom = 0. return endif to
* * use Dalgarno estimate for highly ionized species, number reset with * set charge transfer command if( ipIon.gt.4 ) then *change 96 nov 25, added this option, default is 1.92e-9 HCTRecom = 1.92e-9 * ipIon return endif
This change will have a significant effect on high levels of ionization for models ionized by an X-Ray continuum.
Collision strength for [Ar VI] 4.53 microns was not properly updated. The Iron Project has produced improved collision strengths for many infrared transitions, and a major effort was made to update the code to this work. The collision strength for the [Ar VI] 4.53 micron line was calculated by Saraph and Storey (A&AS 115, 151), as noted in the code, but the collision strength was mistakenly not updated. To correct this change line 267 of arcol.f from
call PutCS( .8 , Ar06453 )
call PutCS( 6.33 , Ar06453 )
Collision strength for [Fe VII] 6087 was too small by a ratio of statistical weights. To correct, add the following after line 545 of routine fecol:
cs = 2.2684 - tused* (0.60893 - 0.052004*tused) *change 97 jan 27, should have been total collision strength of multiplet *thanks for Tino Oliva for finding this problem cs = cs * 2.96 c6087 = popexc(cs ,21.,5.,0.94,2.32e4,fe(7))*3.27E-12 1 * 0.577/0.94
Continuum binning was too coarse near the 4 Ryd He ionization edge. To correct this, change the upper and lower limits on successive calls to routine fill from 4.0 to 4.05 Ryd. Edit routine dopoint, going to lines 395 and 396. Change these two to
*change 97 jan 27, increase resolution between 4 and 4.05 ryd *thanks to Tino Oliva for finding this problem call fill( 3.97 , 4.05 , 0.002, n0 , ipnt ) call fill( 4.05 , 50.0 , 0.008, n0 , ipnt )
Some H- heating cooling derivatives were bad, inhibiting convergence of some models. This did not affect results but could cause convergence problems in low ionization high density clouds where H- heating and cooling are important. Two changes: change line 349 of routine coolr from
dcool = dcool + hmfb * teinv
dcool = dcool + 2.*hmfb * teinv
and in routine SumHeat change line 138 from
dHTotDT = dHTotDT + heating(1,16) * oldfac
dHTotDT = dHTotDT + heating(1,16) * 0.7/te
zones near illuminated face of cloud became too large. The only effect of this is to introduce some noise in predicted line intensities, at the five percent level. This only affected the highest ionization lines, which are formed near the illuminated face. To fix, go to line 672 of routine nextdr and change
if( nzone.lt.11 ) then drmax = 4. * drad else
if( nzone.lt.11 ) then if( nzone.lt.5 ) then drmax = 4. * drad else drmax = 1.3 * drad endif else
Arnaud and Rothenflug collisional ionization rates are discontinuous. The problem was present in all versions of C90, and in any work that uses the Arnaud and Rothenflug representation of the E1 function. The solution is to delete the distributed version of colfit and replace it with the corrected version, present in Dima Verner's Atomic Data for Astrophysics page, or directly by clicking here. The file called colfit.f also contains the routines collis.f, exaife.f, and fcl2.f. All of these original routines must be deleted when the new colfit is downloaded. This problem had little effect on results but could cause temperature failures in some cases.
Acknowledgments: Many thanks to Mike Brotherton, Simon Casassus, John Houck, Anthony Leonard, Valentia Luridiana, Jon Slavin, and Henrik Spoon for finding problems with version 90.02.