123 for( nelem=ipISO; nelem<
LIMELM; ++nelem )
125 HOpacRatSav[ipISO][nelem] = (
double *)
MALLOC(
sizeof(
double)*(
unsigned)(
iso_sp[ipISO][nelem].numLevels_max+1) );
127 hnsav[ipISO][nelem] = (
double *)
MALLOC(
sizeof(
double)*(
unsigned)(
iso_sp[ipISO][nelem].numLevels_max+1) );
136 for(nelem=0; nelem<
LIMELM; ++nelem )
140 (
double*)
MALLOC(
sizeof(
double)*(
unsigned)(nelem+2) );
142 (
double*)
MALLOC(
sizeof(
double)*(
unsigned)(nelem+2) );
145 for( ion=0; ion<nelem+2; ++ion )
157 fprintf(
ioQQQ,
" IterStart called.\n" );
165 rt.lgMaserCapHit =
false;
167 rt.lgMaserSetDR =
false;
174 for(nelem=0; nelem<
HS_NZ; ++nelem )
176 atmdat.lgHCaseBOK[0][nelem] =
true;
177 atmdat.lgHCaseBOK[1][nelem] =
true;
181 for(
long ipSpecies=0; ipSpecies<
nSpecies; ++ipSpecies )
188 strncpy(
StopCalc.chReasonStop,
"reason not specified.",
193 he.frac_he0dest_23S = 0.;
194 he.frac_he0dest_23S_photo = 0.;
199 timesc.time_H2_Dest_longest = 0.;
200 timesc.time_H2_Form_longest = 0.;
202 timesc.time_H2_Dest_here = -1.;
203 timesc.time_H2_Form_here = 0.;
205 timesc.BigCOMoleForm = 0.;
208 hydro.HCollIonMax = 0.;
214 hydro.lgHiPop2 =
false;
245 dense.xMassTotal = 0.;
247 for( nelem=0; nelem <
LIMELM; nelem++ )
250 if(
dense.lgElmtOn[nelem] )
254 for( ion=0; ion < (nelem + 2); ion++ )
259 for( ion=0; ion <
LIMELM; ion++ )
273 for( nelem=ipISO; nelem <
LIMELM; nelem++ )
275 if(
dense.lgElmtOn[nelem] )
277 for( n=0; n <
iso_sp[ipISO][nelem].numLevels_max; n++ )
280 hnsav[ipISO][nelem][n] =
iso_sp[ipISO][nelem].st[n].Pop();
283 for( vector<two_photon>::iterator tnu =
iso_sp[ipISO][nelem].TwoNu.begin(); tnu !=
iso_sp[ipISO][nelem].TwoNu.end(); ++tnu )
284 tnu->induc_dn_max = 0.;
288 rfield.ipEnergyBremsThin = 0;
289 rfield.EnergyBremsThin = 0.;
321 for( ion=0; ion<nelem+2; ++ion )
329 for( ion2=0; ion2<nelem+2; ++ion2 )
362 radius.dr_max_last_iter = 0.;
375 rfield.extin_mag_B_point = 0.;
376 rfield.extin_mag_V_point = 0.;
377 rfield.extin_mag_B_extended = 0.;
378 rfield.extin_mag_V_extended = 0.;
381 for( i=0; i <
rfield.nupper+1; i++ )
386 rfield.ConEmitReflec[0][i] = 0.;
387 rfield.ConEmitOut[0][i] = 0.;
388 rfield.ConInterOut[i] = 0.;
393 rfield.outlin_noplot[i] = 0.;
394 rfield.ConOTS_local_OTS_rate[i] = 0.;
395 rfield.ConOTS_local_photons[i] = 0.;
396 rfield.DiffuseEscape[i] = 0.;
399 opac.opacity_abs_savzon1[i] =
opac.opacity_abs[i];
400 opac.opacity_sct_savzon1[i] =
opac.opacity_sct[i];
431 for( i=0; i <
NCOLD; i++ )
438 colden.coldenH2_ov_vel = 0.;
446 (*diatom)->ortho_colden = 0.;
447 (*diatom)->para_colden = 0.;
452 mole.species[i].column = 0.;
454 for( i=0; i < 5; i++ )
462 for( i=0; i < 3; i++ )
473 for( i=0; i < 4; i++ )
488 for(
unsigned i = 0; i <
mole.species.size(); ++i )
490 if(
mole.species[i].levels != NULL )
519 for( ion=0; ion<nelem+1; ++ion )
529 for( ion=0; ion<nelem+1; ++ion )
531 ionbal.CompRecoilHeatRateSave[nelem][ion] =
ionbal.CompRecoilHeatRate[nelem][ion];
532 ionbal.CompRecoilIonRateSave[nelem][ion] =
ionbal.CompRecoilIonRate[nelem][ion];
538 conv.nTotalFailures = 0;
558 hmi.HeatH2DexcMax = 0.;
559 hmi.CoolH2DexcMax = 0.;
561 hmi.h2line_cool_frac = 0.;
586 LineSv[i].SumLine[0] = 0.;
587 LineSv[i].SumLine[1] = 0.;
608 hmi.h2line_cool_frac = 0.;
621 fprintf(
ioQQQ,
"PROBLEM could not find the normalisation line.\n");
623 "IterStart could not find a line with a wavelength of %f and a label of %s\n",
625 fprintf(
ioQQQ,
"Please check the emission line output to find the correct line identification.\n");
626 fprintf(
ioQQQ,
"Sorry.\n");
635 LineSave.ipNormWavL =
cdLine(
"TOTL" , 4861.36f , &fhe1nx , &ratio );
639 fprintf(
ioQQQ,
"PROBLEM could not find the default normalisation line.\n");
641 "IterStart could not find a line with a wavelength of 4861 and a label of TOTL\n" );
642 fprintf(
ioQQQ,
"Please check the emission line output to find the correct line identification.\n");
643 fprintf(
ioQQQ,
"Sorry.\n");
654 for(
long int nStopLine=0; nStopLine <
StopCalc.nstpl; nStopLine++ )
656 double relint, absint ;
661 StopCalc.StopLineWl1[nStopLine], &relint, &absint );
663 if(
StopCalc.ipStopLin1[nStopLine]<0 )
666 " IterStart could not find first line in STOP LINE command, number %ld with label *%s* and wl %.1f\n",
675 StopCalc.StopLineWl2[nStopLine], &relint, &absint );
677 if(
StopCalc.ipStopLin2[nStopLine] < 0 )
680 " IterStart could not find second line in STOP LINE command, line number %ld with label *%s* and wl %.1f\n",
690 " stop line 1 is number %5ld wavelength is %f label is %4.4s\n",
695 " stop line 2 is number %5ld wavelength is %f label is %4.4s\n",
704 if(
prt.lgPrtLastIt )
725 fprintf(
ioQQQ,
" IterStart does not change mid-zone optical depth since CASE B\n" );
744 if( ratio >
hydro.FracInd )
754 if( ratio >
hydro.fbul )
762 fprintf(
ioQQQ,
" IterStart returns.\n" );
789 if(
StopCalc.TeFloor > 0. && !
thermal.lgTemperatureConstantCommandParsed )
791 thermal.lgTemperatureConstant =
false;
804 (*diatom)->H2_Reset();
809 for( ion=0; ion<nelem+2; ++ion )
817 for( ion2=0; ion2<nelem+2; ++ion2 )
829 dense.updateXMolecules();
842 for( i=0; i <
NCOLD; i++ )
847 colden.coldenH2_ov_vel = 0.;
852 mole.species[i].column = 0.;
854 mole.species[i].xFracLim = 0.;
855 mole.species[i].nAtomLim = -1;
863 rfield.extin_mag_B_point = 0.;
864 rfield.extin_mag_V_point = 0.;
865 rfield.extin_mag_B_extended = 0.;
866 rfield.extin_mag_V_extended = 0.;
886 rfield.ipEnergyBremsThin = 0;
887 rfield.EnergyBremsThin = 0.;
897 radius.lgDrMinUsed =
false;
908 if(
trace.nTrConvg==0 )
910 trace.lgTrace =
true;
917 fprintf(
ioQQQ,
" IterRestart called.\n" );
921 trace.lgTrace =
false;
927 for( ion=0; ion<nelem+1; ++ion )
935 for( nelem=0; nelem<
LIMELM; ++nelem)
937 for( ion=0; ion<nelem+1; ++ion )
939 ionbal.CompRecoilHeatRate[nelem][ion] =
ionbal.CompRecoilHeatRateSave[nelem][ion];
940 ionbal.CompRecoilIonRate[nelem][ion] =
ionbal.CompRecoilIonRateSave[nelem][ion];
944 wind.lgVelPos =
true;
988 for(
long nelem=ipISO; nelem <
LIMELM; ++nelem )
990 iso_sp[ipISO][nelem].Reset();
996 for( nelem=ipISO; nelem <
LIMELM; nelem++ )
998 if(
dense.lgElmtOn[nelem] )
1000 for( n=
ipH1s; n <
iso_sp[ipISO][nelem].numLevels_max; n++ )
1003 iso_sp[ipISO][nelem].st[n].Pop() =
hnsav[ipISO][nelem][n];
1011 fprintf(
ioQQQ,
" EDEN set to%12.4e by IterRestart.\n",
1017 for( ion=0; ion < (nelem + 2); ion++ )
1021 for( i=0; i <
LIMELM; i++ )
1032 rfield.resetCoarseTransCoef();
1035 for( i=0; i <
rfield.nupper; i++ )
1038 rfield.flux_beam_const[i] =
rfield.flux_beam_const_save[i];
1041 rfield.flux_beam_time[i] =
rfield.flux_time_beam_save[i]*
1042 rfield.time_continuum_scale;
1050 if( fac_old > log10(DBL_MAX) )
1054 else if( fac_new > log10(DBL_MAX) )
1058 else if( fac_old > 1.e-5 )
1060 slope_ratio = (exp(fac_old) - 1.)/(exp(fac_new) - 1.);
1064 slope_ratio = (fac_old*(1. + fac_old/2.))/(fac_new*(1. + fac_new/2.));
1072 rfield.flux_isotropic[i] =
rfield.flux_isotropic_save[i];
1076 rfield.flux_beam_const[i];
1085 rfield.ConInterOut[i] = 0.;
1086 rfield.OccNumbBremsCont[i] = 0.;
1087 rfield.OccNumbDiffCont[i] = 0.;
1088 rfield.OccNumbContEmitOut[i] = 0.;
1089 rfield.outlin[0][i] = 0.;
1090 rfield.outlin_noplot[i] = 0.;
1091 rfield.ConOTS_local_OTS_rate[i] = 0.;
1092 rfield.ConOTS_local_photons[i] = 0.;
1093 rfield.DiffuseEscape[i] = 0.;
1095 opac.opacity_abs[i] =
opac.opacity_abs_savzon1[i];
1096 opac.OldOpacSave[i] =
opac.opacity_abs_savzon1[i];
1097 opac.opacity_sct[i] =
opac.opacity_sct_savzon1[i];
1120 rfield.reflin[0][i] = 0.;
1121 rfield.ConEmitReflec[0][i] = 0.;
1122 rfield.ConEmitOut[0][i] = 0.;
1123 rfield.ConRefIncid[0][i] = 0.;
1135 opac.E2TauAbsOut[i] = 1.;
1148 fprintf(
ioQQQ,
"\f\n Start Iteration Number %ld %75.75s\n\n\n",
1185 double cumulative_factor = (
dynamics.timestep/
1188 for(
long n=0; n<
LineSave.nsum; ++n )
1190 for(
long nEmType=0; nEmType<2; ++nEmType )
1196 for(
long n=0; n<
rfield.nflux; ++n)
1203 rfield.flux_total_incident[1][n] += (
realnum)
rfield.flux_total_incident[0][n]*cumulative_factor;
1212 for(
long i=0; i<
struc.nzonePreviousIteration; ++i )
1221 for(
long i=0; i <
rfield.nflux; i++ )
1224 enum{DEBUG_LOC=
false};
1227 fprintf(
ioQQQ,
"i=%li opac %.2e \n", i,
1228 (
double)
opac.opacity_abs[i]*
radius.drad_x_fillfac/2.*(
double)
geometry.DirectionalCosin );