cloudy  trunk
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ion_recomb.cpp
Go to the documentation of this file.
1 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and
2  * others. For conditions of distribution and use see copyright notice in license.txt */
3 /*ion_recomb generate recombination coefficients for any species */
4 /*ion_recombAGN generate recombination coefficients for AGN table */
5 #include "cddefines.h"
6 #include "phycon.h"
7 #include "heavy.h"
8 #include "hmi.h"
9 #include "grainvar.h"
10 #include "dense.h"
11 #include "conv.h"
12 #include "thermal.h"
13 #include "iso.h"
14 #include "abund.h"
15 #include "punch.h"
16 #include "elementnames.h"
17 #include "atmdat.h"
18 #include "ionbal.h"
19 
20 /*ion_recomb generate recombination coefficients for any species */
22  /* this is debug flag */
23  bool lgPrintIt,
24  const double *dicoef,
25  const double *dite,
26  const double ditcrt[],
27  const double aa[],
28  const double bb[],
29  const double cc[],
30  const double dd[],
31  const double ff[],
32  /* nelem is the atomic number on the C scale, 0 for H */
33  long int nelem/*,
34  double tlow[]*/)
35 {
36 #define DICOEF(I_,J_) (*(dicoef+(I_)*(nelem+1)+(J_)))
37 #define DITE(I_,J_) (*(dite+(I_)*(nelem+1)+(J_)))
38  long int i,
39  ion,
40  limit;
41  double
42  fac2,
43  factor,
44  DielRecomRateCoef_HiT[LIMELM] ,
45  DielRecomRateCoef_LowT[LIMELM] ,
46  ChargeTransfer[LIMELM] ,
47  t4m1,
48  tefac;
49 
50  /* these are used for adding noise to rec coef */
51  static double RecNoise[LIMELM][LIMELM];
52  static bool lgNoiseNeedEval=true;
53 
54  DEBUG_ENTRY( "ion_recomb(false,)" );
55 
56  /* set up ionization balance matrix, C(I,1)=destruction, 2=creation
57  * heating rates saved in array B(I) in same scratch block
58  * factor is for Aldrovandi+Pequignot fit, FAC2 is for Nuss+Storey
59  * fit for dielectronic recombination
60  * GrnIonRec is rate ions recombine on grain surface, normally zero;
61  * set in hmole, already has factor of hydrogen density
62  * */
63 
64  /* routine only used for Li on up */
65  ASSERT( nelem < LIMELM);
66  ASSERT( nelem > 1 );
67 
68  /* check that range of ionization is correct */
69  ASSERT( dense.IonLow[nelem] >= 0 );
70  ASSERT( dense.IonLow[nelem] <= nelem+1 );
71 
73  t4m1 = 1e4/phycon.te;
74  fac2 = 1e-14*phycon.sqrte;
75 
76  /* option to put noise into rec coefficient -
77  * one time initialization of noise - this is set with command
78  * set dielectronic recombination kludge noise */
79  if( lgNoiseNeedEval )
80  {
81  int n;
82  if( ionbal.guess_noise !=0. )
83  {
84  for( n=ipHYDROGEN; n<LIMELM; ++n )
85  {
86  for( ion=0; ion<=n; ++ion )
87  {
88  /* log normal noise with dispersion entered on command line */
89  /* NB the seed for rand was set when the command was parsed */
90  RecNoise[n][ion] = pow(10., RandGauss( 0. , ionbal.guess_noise ) );
91  }
92  }
93  }
94  else
95  {
96  for( n=ipHYDROGEN; n<LIMELM; ++n )
97  {
98  for( ion=0; ion<=n; ++ion )
99  {
100  RecNoise[n][ion] = 1.;
101  }
102  }
103  }
104  lgNoiseNeedEval = false;
105  }
106 
107  /* this routine only does simple two-level species,
108  * loop over ions will be <= limit, IonHigh is -1 since very
109  * highest stage of ionization is not recombined into.
110  * for Li, will do only atom, since ions are H and He like,
111  * so limit is zero */
112  limit = MIN2(nelem-NISO,dense.IonHigh[nelem]-1);
113  ASSERT( limit >= -1 );
114 
115  /* zero-out loop comes before main loop since there are off-diagonal
116  * elements in the main ionization loop, due to multi-electron processes */
117  /* >>chng 00 dec 07, limit changed to identical to ion_solver */
118  for( ion=0; ion <= limit; ion++ )
119  {
120  ionbal.RateRecomTot[nelem][ion] = 0.;
121  ChargeTransfer[ion] = 0.;
122  DielRecomRateCoef_LowT[ion] = 0.;
123  DielRecomRateCoef_HiT[ion] = 0.;
124  ionbal.RR_rate_coef_used[nelem][ion] = 0.;
125  ionbal.DR_rate_coef_used[nelem][ion] = 0.;
126  }
127  for( ion=limit+1; ion < LIMELM; ion++ )
128  {
129  /* >>chng 01 dec 18, do not set this to FLT_MAX since it clobbers what
130  * had been set in h-like and he-like routines - that would only affect
131  * the printout */
132  ChargeTransfer[ion] = -FLT_MAX;
133  DielRecomRateCoef_LowT[ion] = -FLT_MAX;
134  DielRecomRateCoef_HiT[ion] = -FLT_MAX;
135  }
136 
137  DielRecomRateCoef_HiT[nelem] = 0.;
138  DielRecomRateCoef_HiT[nelem-1] = 0.;
139 
140  DielRecomRateCoef_LowT[nelem] = 0.;
141  DielRecomRateCoef_LowT[nelem-1] = 0.;
142 
143  /* these are counted elsewhere and must not be added here */
144  Heavy.xLyaHeavy[nelem][nelem] = 0.;
145  Heavy.xLyaHeavy[nelem][nelem-1] = 0.;
146 
147  /* IonLow is 0 for the atom, limit chosen to NOT include iso sequences */
148  for( ion=dense.IonLow[nelem]; ion <= limit; ion++ )
149  {
150  /* number of bound electrons of the ion after recombination,
151  * for an atom (ion=0) this is
152  * equal to nelem+1, the element on the physical scale, since nelem is
153  * on the C scale, being 5 for carbon */
154  long n_bnd_elec_after_recomb = nelem+1 - ion;
155 
156  /* >>chng 02 nov 06 add charge transfer here rather than in ion_solver */
157  /* charge transfer recombination of this species by ionizing hydrogen and helium */
158  ChargeTransfer[ion] =
159  /* He0 + ion charge transfer recombination */
160  atmdat.HeCharExcRecTo[nelem][ion]*
161  /* following product is density [cm-3] of ground state of He0 */
163  /* H0 + ion charge transfer recombination */
164  atmdat.HCharExcRecTo[nelem][ion]*
165  /* following product is density [cm-3] of ground state of H0 */
167 
168  /*>>chng 04 feb 20, add this, had always been in for destruction for H- */
169  /* charge transfer recombination of first ion to neutral, by reaction with H-
170  * the ion==0 is right, the first array element is the */
171  if( ion==0 && nelem>ipHELIUM && atmdat.lgCTOn )
172  ChargeTransfer[ion] += hmi.hmin_ct_firstions * hmi.Hmolec[ipMHm];
173 
174  /* >>chng 06 feb 01, add option to use Badnell RR data rather than Verner */
176  {
177  ionbal.RR_rate_coef_used[nelem][ion] = ionbal.RR_Badnell_rate_coef[nelem][ion];
178  }
179  else
180  {
181  ionbal.RR_rate_coef_used[nelem][ion] = ionbal.RR_Verner_rate_coef[nelem][ion];
182  }
183 
184  /* Burgess or high-T dielectronic recombination */
185  DielRecomRateCoef_HiT[ion] = 0.;
186  /* >>chng 03 oct 29, do fe here rather than after loop */
187  if( nelem==ipIRON )
188  {
189  /* implements dn > 0 DR from Arnaud & Raymond 1992 */
190  DielRecomRateCoef_HiT[ion] = atmdat_dielrec_fe(ion+1,phycon.te);
191  }
192  else if( phycon.te > (ditcrt[ion]*0.1) )
193  {
194  DielRecomRateCoef_HiT[ion] = ionbal.DielSupprs[0][ion]/phycon.te32*
195  DICOEF(0,ion)*exp(-DITE(0,ion)/phycon.te)*
196  (1. + DICOEF(1,ion)*
197  sexp(DITE(1,ion)/phycon.te));
198  }
199 
200  /* begin dn = 0 dielectronic recombination
201  * do not include it for rec from
202  * a closed shell, n_bnd_elec_after_recomb-1 is number of bound electrons in parent ion */
203  DielRecomRateCoef_LowT[ion] = 0.;
204  if( ((n_bnd_elec_after_recomb-1) != 2) &&
205  ((n_bnd_elec_after_recomb-1) != 10) &&
206  ((n_bnd_elec_after_recomb-1) != 18) )
207  {
208  tefac = ff[ion]*t4m1;
209  /* do not do iron here since all dn = 0 DR either Badnell or a guess */
210  if( ff[ion] != 0. && nelem!=ipIRON )
211  {
212  /* >>chng 06 feb 14, O+3 ff[ion] is very negative, as a result exp goes to +inf
213  * at very low temperatures. This is a error in the Nussbaumer & Storey fits to DR.
214  * do not use them is tefac = ff[ion] / t4 is very negative */
215  if( tefac > -5. )
216  {
217  factor = (((aa[ion]*t4m1+bb[ion])*t4m1+cc[ion])*t4m1+dd[ion])* sexp(tefac);
218  DielRecomRateCoef_LowT[ion] = ionbal.DielSupprs[1][ion]*fac2*
219  MAX2(0.,factor );
220  }
221  else
222  {
223  DielRecomRateCoef_LowT[ion] = 0.;
224  }
225  }
226  else if( ionbal.lg_guess_coef )
227  {
228  /* first guesses are those based on Nussbaumer & Storey and are from
229  * >>refer ion DR Ali, B., Blum, R. D., Bumgardner, T. E., Cranmer, S. R.,
230  * >>refercon Ferland, G. J., Haefner, R. I., & Tiede, G. P. 1991, PASP, 103, 1182*/
231  if( (ff[ion] == 0.) && (ion <= 3) )
232  {
233  /* these are guessed dielectronic rec rates, taken from
234  * >>refer all dielectronic Ali, B., Blum, R. D., Bumgardner, T. E., Cranmer, S. R., Ferland, G. J.,
235  * >>refercon Haefner, R. I., & Tiede, G. P. 1991, PASP, 103, 1182 */
236  static double cludge[4]={3e-13,3e-12,1.5e-11,2.5e-11};
237 
238  /* >>chng 01 jun 30, make GuessDiel array */
239  DielRecomRateCoef_LowT[ion] = ionbal.DielSupprs[1][ion]*cludge[ion]*
240  /* this is optional scale factor on set dielectronic rec command
241  * >>chng 06 feb 07, add Boltzmann factor to induce behavior
242  * like in Nussbaumer & Storey */
243  ionbal.GuessDiel[ion] * sexp( 1e3 / phycon.te );
244  }
245  /* this implements the low-T kludge dielectronic command */
246  else
247  {
248  /* assume the recombination rate is the coefficient scanned off the steve option, times the charge
249  * before recombination raised to a power */
250  double fitcoef[3][3] =
251  {
252  /* L- shell, Z=17-23 */
253  {-52.5073,+5.19385,-0.126099} ,
254  /* M-shell (Z=9-15) */
255  {-10.9679,1.66397,-0.0605965} ,
256  /* N shell */
257  {-3.95599,1.61884,-0.194540}
258  };
259 
260  /* these implement the guesses in
261  * Kraemer, S.B., Ferland, G.J., & Gabel, J.R. 2004, ApJ, 604, 556 */
262  if( nelem==ipIRON )
263  {
264  int nshell;
265  if( (n_bnd_elec_after_recomb>=4) && (n_bnd_elec_after_recomb<=11) )
266  nshell = 0;
267  else if( (n_bnd_elec_after_recomb>=12) && (n_bnd_elec_after_recomb<=19 ))
268  nshell = 1;
269  else
270  nshell = 2;
271  /* n_bnd_elec_after_recomb is the number of bound electrons */
272  DielRecomRateCoef_LowT[ion] = fitcoef[nshell][0] +
273  fitcoef[nshell][1]*(ion+1) +
274  fitcoef[nshell][2]*POW2(ion+1.);
275  DielRecomRateCoef_LowT[ion] = 1e-10*pow(10.,DielRecomRateCoef_LowT[ion]);
276 
277  }
278  else
279  {
280  /* this is guess for all other elements presented in
281  * >>refer all DR Kraemer, S.B., Ferland, G.J., & Gabel, J.R. 2004, ApJ, 604, 556 */
282  DielRecomRateCoef_LowT[ion] = 3e-12*pow(10.,(double)(ion+1)*0.1);
283  }
284  /* >>chng 06 feb 02, add option to use mean of Badnell DR in place
285  * of these hacks */
287  DielRecomRateCoef_LowT[ion] = ionbal.DR_Badnell_rate_coef_mean_ion[ion];
288  }
289  /* include optional noise here
290  * >>chng 06 feb 07, move noise down to here so that use for both
291  * guesses of DR rates */
292  DielRecomRateCoef_LowT[ion] *= RecNoise[nelem][ion];
293  }
294  }
295  /* >>chng 05 dec 19, add option to use Badnell numbers */
296  /* this is total old DR rates - may not use it */
297  ionbal.DR_old_rate_coef[nelem][ion] = DielRecomRateCoef_HiT[ion] + DielRecomRateCoef_LowT[ion];
298 
299  /* set total DR rate - either Badnell if it exists or on with set dielectronic recombination badnell */
301  {
302  ionbal.DR_rate_coef_used[nelem][ion] = ionbal.DR_Badnell_rate_coef[nelem][ion];
303  }
304  else
305  {
306  ionbal.DR_rate_coef_used[nelem][ion] = ionbal.DR_old_rate_coef[nelem][ion];
307  }
308 
309  /* sum of recombination rates [units s-1] for radiative, three body, charge transfer */
310  ionbal.RateRecomTot[nelem][ion] =
311  dense.eden* (
312  ionbal.RR_rate_coef_used[nelem][ion] +
313  ionbal.DR_rate_coef_used[nelem][ion] +
314  ionbal.CotaRate[ion] ) +
315  ChargeTransfer[ion];
316 
317  /* >>chng 01 jun 30, FRAC_LINE was 0.1, not 1, did not include anything except
318  * radiative recombination, the radrec term */
319 # define FRAC_LINE 1.
320  /* was 0.1 */
321  /*Heavy.xLyaHeavy[nelem][ion] = (realnum)(dense.eden*radrec*FRAC_LINE );*/
322  Heavy.xLyaHeavy[nelem][ion] = (realnum)(dense.eden*
323  (ionbal.RR_rate_coef_used[nelem][ion]+DielRecomRateCoef_LowT[ion]+DielRecomRateCoef_HiT[ion])*FRAC_LINE );
324  }
325 
326  /* option to punch rec coefficients */
327  if( punch.lgioRecom || lgPrintIt )
328  {
329  /* >>chng 04 feb 22, make option to print ions for single element */
330  FILE *ioOut;
331  if( lgPrintIt )
332  ioOut = ioQQQ;
333  else
334  ioOut = punch.ioRecom;
335 
336  /* print name of element */
337  fprintf( ioOut,
338  " %s recombination coefficients fnzone:%.2f \tte\t%.4e\tne\t%.4e\n",
340 
341  /*limit = MIN2(11,dense.IonHigh[nelem]);*/
342  /* >>chng 05 sep 24, just print one long line - need info */
343  limit = dense.IonHigh[nelem];
344  for( i=0; i < limit; i++ )
345  {
346  fprintf( ioOut, "%10.2e", ionbal.RR_rate_coef_used[nelem][i] );
347  }
348  fprintf( ioOut, " radiative used vs Z\n" );
349 
350  for( i=0; i < limit; i++ )
351  {
352  fprintf( ioOut, "%10.2e", ionbal.RR_Verner_rate_coef[nelem][i] );
353  }
354  fprintf( ioOut, " old Verner vs Z\n" );
355 
356  for( i=0; i < limit; i++ )
357  {
358  fprintf( ioOut, "%10.2e", ionbal.RR_Badnell_rate_coef[nelem][i] );
359  }
360  fprintf( ioOut, " new Badnell vs Z\n" );
361 
362  for( i=0; i < limit; i++ )
363  {
364  /* >>chng 06 jan 19, from div by eden to div by H0 - want units of cm3 s-1 but
365  * no single collider does this so not possible to get rate coefficient easily
366  * H0 is more appropriate than electron density */
367  fprintf( ioOut, "%10.2e", ChargeTransfer[i]/SDIV(dense.xIonDense[ipHYDROGEN][0]) );
368  }
369  fprintf( ioOut, " CT/n(H0)\n" );
370 
371  for( i=0; i < limit; i++ )
372  {
373  fprintf( ioOut, "%10.2e", ionbal.CotaRate[ion] );
374  }
375  fprintf( ioOut, " 3body vs Z /ne\n" );
376 
377  /* note different upper limit - this routine does grain rec for all ions */
378  for( i=0; i < dense.IonHigh[nelem]; i++ )
379  {
380  fprintf( ioOut, "%10.2e", gv.GrainChTrRate[nelem][i+1][i]/dense.eden );
381  }
382  fprintf( ioOut, " Grain vs Z /ne\n" );
383 
384  for( i=0; i < limit; i++ )
385  {
386  fprintf( ioOut, "%10.2e", DielRecomRateCoef_HiT[i] );
387  }
388  fprintf( ioOut, " Burgess vs Z\n" );
389 
390  for( i=0; i < limit; i++ )
391  {
392  fprintf( ioOut, "%10.2e", ionbal.DR_old_rate_coef[nelem][i] );
393  }
394  fprintf( ioOut, " old Nussbaumer Storey DR vs Z\n" );
395 
396  for( i=0; i < limit; i++ )
397  {
398  fprintf( ioOut, "%10.2e", ionbal.DR_Badnell_rate_coef[nelem][i] );
399  }
400  fprintf( ioOut, " new Badnell DR vs Z\n" );
401 
402  for( i=0; i < limit; i++ )
403  {
404  fprintf( ioOut, "%10.2e", DielRecomRateCoef_LowT[i] );
405  }
406  fprintf( ioOut, " low T DR used vs Z\n" );
407 
408  /* total recombination rate, with density included - this goes into the matrix */
409  for( i=0; i < limit; i++ )
410  {
411  fprintf( ioOut, "%10.2e", ionbal.RateRecomTot[nelem][i] );
412  }
413  fprintf( ioOut,
414  " total rec rate (with density) for %s\n",
415  elementnames.chElementSym[nelem] );
416  for( i=0; i < limit; i++ )
417  {
418  fprintf( ioOut, "%10.2e", ionbal.RateRecomTot[nelem][i]/dense.eden );
419  }
420  fprintf( ioOut,
421  " total rec rate / ne for %s\n\n",
422  elementnames.chElementSym[nelem] );
423 
424  /* spill over to next line for many stages of ionization */
425  if( dense.IonHigh[nelem] > 11 )
426  {
427  limit = MIN2(29,dense.IonHigh[nelem]);
428  fprintf( ioOut, " R " );
429  for( i=11; i < limit; i++ )
430  {
431  fprintf( ioOut, "%10.2e", dense.eden*ionbal.CotaRate[ion] );
432  }
433  fprintf( ioOut, "\n" );
434 
435  fprintf( ioOut, " B " );
436  for( i=11; i < limit; i++ )
437  {
438  fprintf( ioOut, "%10.2e", DielRecomRateCoef_HiT[i] );
439  }
440  fprintf( ioOut, "\n" );
441 
442  fprintf( ioOut, " NS" );
443  for( i=11; i < limit; i++ )
444  {
445  fprintf( ioOut, "%10.2e", DielRecomRateCoef_LowT[i] );
446  }
447  fprintf( ioOut, "\n" );
448 
449  fprintf( ioOut, " " );
450  for( i=11; i < limit; i++ )
451  {
452  fprintf( ioOut, "%10.2e", ionbal.RateRecomTot[nelem][i] );
453  }
454  fprintf( ioOut, "\n\n" );
455  }
456  }
457 
458  /* >>chng 02 nov 09, from -2 to -NISO */
459  /*limit = MIN2(nelem-2,dense.IonHigh[nelem]-1);*/
460  limit = MIN2(nelem-NISO,dense.IonHigh[nelem]-1);
461  for( i=dense.IonLow[nelem]; i <= limit; i++ )
462  {
463  ASSERT( Heavy.xLyaHeavy[nelem][i] > 0. );
464  ASSERT( ionbal.RateRecomTot[nelem][i] > 0. );
465  }
466  return;
467 # undef DITE
468 # undef DICOEF
469 }
470 
471 /*ion_recombAGN generate recombination coefficients for AGN table */
472 void ion_recombAGN( FILE * io )
473 {
474 # define N1LIM 3
475 # define N2LIM 4
476  double te1[N1LIM]={ 5000., 10000., 20000.};
477  double te2[N2LIM]={ 20000.,50000.,100000.,1e6};
478  /* this is boundary between two tables */
479  double BreakEnergy = 100./13.0;
480  long int nelem, ion , i;
481  /* this will hold element symbol + ionization */
482  char chString[100],
483  chOutput[100];
484  /* save temp here */
485  double TempSave = phycon.te;
486  /* save ne here */
487  double EdenSave = dense.eden;
488 
489  DEBUG_ENTRY( "ion_recomb(false,)" );
490 
491  dense.eden = 1.;
492  /*atmdat_readin();*/
493 
494  /* first put header on file */
495  fprintf(io,"X+i\\Te");
496  for( i=0; i<N1LIM; ++i )
497  {
498  phycon.te = te1[i];
499  fprintf(io,"\t%.0f K",phycon.te);
500  }
501  fprintf(io,"\n");
502 
503  /* now do loop over temp, but add elements */
504  for( nelem=ipLITHIUM; nelem<LIMELM; ++nelem )
505  {
506  /* this list of elements included in the AGN tables is defined in zeroabun.c */
507  if( abund.lgAGN[nelem] )
508  {
509  for( ion=0; ion<=nelem; ++ion )
510  {
511  ASSERT( Heavy.Valence_IP_Ryd[nelem][ion] > 0.05 );
512 
513  if( Heavy.Valence_IP_Ryd[nelem][ion] > BreakEnergy )
514  break;
515 
516  /* print chemical symbol */
517  sprintf(chOutput,"%s",
518  elementnames.chElementSym[nelem]);
519  /* some elements have only one letter - this avoids leaving a space */
520  if( chOutput[1]==' ' )
521  chOutput[1] = chOutput[2];
522  /* now ionization stage */
523  if( ion==0 )
524  {
525  sprintf(chString,"0 ");
526  }
527  else if( ion==1 )
528  {
529  sprintf(chString,"+ ");
530  }
531  else
532  {
533  sprintf(chString,"+%li ",ion);
534  }
535  strcat( chOutput , chString );
536  fprintf(io,"%5s",chOutput );
537 
538  for( i=0; i<N1LIM; ++i )
539  {
540  TempChange(te1[i] , false);
541  dense.IonLow[nelem] = 0;
542  dense.IonHigh[nelem] = nelem+1;
543  if( ConvBase(0) )
544  fprintf(ioQQQ,"PROBLEM ConvBase returned error.\n");
545  fprintf(io,"\t%.2e",ionbal.RateRecomTot[nelem][ion]);
546  }
547  fprintf(io,"\n");
548  }
549  fprintf(io,"\n");
550  }
551  }
552 
553  /* second put header on file */
554  fprintf(io,"X+i\\Te");
555  for( i=0; i<N2LIM; ++i )
556  {
557  TempChange(te2[i] , false);
558  fprintf(io,"\t%.0f K",phycon.te);
559  }
560  fprintf(io,"\n");
561 
562  /* now do same loop over temp, but add elements */
563  for( nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
564  {
565  /* this list of elements included in the AGN tables is defined in zeroabun.c */
566  if( abund.lgAGN[nelem] )
567  {
568  for( ion=0; ion<=nelem; ++ion )
569  {
570  ASSERT( Heavy.Valence_IP_Ryd[nelem][ion] > 0.05 );
571 
572  if( Heavy.Valence_IP_Ryd[nelem][ion] <= BreakEnergy )
573  continue;
574 
575  /* print chemical symbol */
576  fprintf(io,"%s",
577  elementnames.chElementSym[nelem]);
578  /* now ionization stage */
579  if( ion==0 )
580  {
581  fprintf(io,"0 ");
582  }
583  else if( ion==1 )
584  {
585  fprintf(io,"+ ");
586  }
587  else
588  {
589  fprintf(io,"+%li",ion);
590  }
591 
592  for( i=0; i<N2LIM; ++i )
593  {
594  TempChange(te2[i] , false);
595  dense.IonLow[nelem] = 0;
596  dense.IonHigh[nelem] = nelem+1;
597  if( ConvBase(0) )
598  fprintf(ioQQQ,"PROBLEM ConvBase returned error.\n");
599  fprintf(io,"\t%.2e",ionbal.RateRecomTot[nelem][ion]);
600  }
601  fprintf(io,"\n");
602  }
603  fprintf(io,"\n");
604  }
605  }
606 
607  TempChange(TempSave , true);
608  dense.eden = EdenSave;
609  return;
610 }

Generated for cloudy by doxygen 1.8.4