60 static int nCalled = 0;
77 fprintf(
ioQQQ,
" ContCreatePointers called, not evaluating.\n" );
83 fprintf(
ioQQQ,
" ContCreatePointers called first time.\n" );
87 for(
long i=0; i <
rfield.nupper; i++ )
90 strcpy(
rfield.chContLabel[i],
" ");
91 strcpy(
rfield.chLineLabel[i],
" ");
97 for(
long nelem=0; nelem<
LIMELM; ++nelem )
99 if(
dense.lgElmtOn[nelem] )
101 for(
long ion=0; ion<
LIMELM; ++ion )
103 for(
long nshells=0; nshells<7; ++nshells )
105 for(
long j=0; j<3; ++j )
107 opac.ipElement[nelem][ion][nshells][j] = INT_MIN;
130 for(
long nelem=ipISO; nelem < 2; nelem++ )
133 sprintf( chLab,
"%2s%2ld",
elementnames.chElementSym[nelem], nelem+1-ipISO );
137 for(
long ipHi=1; ipHi <
iso_sp[ipISO][nelem].numLevels_max; ipHi++ )
143 for(
long ipLo=0; ipLo < ipHi; ipLo++ )
145 if(
iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().Aul() <=
iso_ctrl.SmallA )
152 if(
iso_sp[ipISO][nelem].trans(ipHi,ipLo).EnergyRyd() <
continuum.filbnd[0] )
156 iso_sp[ipISO][nelem].trans(ipHi,ipLo).ipCont() =
158 iso_sp[ipISO][nelem].fb[ipLo].ipIsoLevNIonCon);
159 iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().ipFine() =
163 if(
iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().ipFine() > 0 )
168 realnum widFine = anuFine -
rfield.fine_anu[
iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().ipFine()-1];
175 ASSERT( fabs(anuCoarse - anuFine) / anuCoarse <
203 fprintf(
ioQQQ,
" insanity heii pointer fix contcreatepointers\n");
294 for(
long i=1; i<
rfield.nupper-1; ++i )
324 for(
long i=0; i <
rfield.nupper; i++ )
341 Heavy.nsShells[nelem][0] = 1;
348 opac.ipElement[nelem][ion][0][1] =
rfield.nupper;
358 Heavy.nsShells[nelem][0] = 1;
365 opac.ipElement[nelem][ion][0][1] =
rfield.nupper;
370 Heavy.nsShells[nelem][1] = 1;
377 opac.ipElement[nelem][ion][0][1] =
rfield.nupper;
397 for(
long nelem=2; nelem<
LIMELM; ++nelem )
401 if(
dense.lgElmtOn[nelem])
404 for(
long j=0; j<=nelem; ++j )
408 for(
long j=0; j<nelem; ++j )
416 for(
long nelem=0; nelem<
LIMELM; ++nelem )
418 if(
dense.lgElmtOn[nelem])
420 for(
long ion=0; ion<nelem+1; ++ion )
468 fprintf(
ioQQQ,
" ContCreatePointers:%ld energy cells used. N(1R):%4ld N(1.8):%4ld N(4Ryd):%4ld N(O3)%4ld N(x-ray):%5ld N(rcoil)%5ld\n",
478 fprintf(
ioQQQ,
" ContCreatePointers: ipEnerGammaRay: %5ld IPPRpari produc%5ld\n",
481 fprintf(
ioQQQ,
" ContCreatePointers: H pointers;" );
482 for(
long i=0; i <= 6; i++ )
486 fprintf(
ioQQQ,
"\n" );
488 fprintf(
ioQQQ,
" ContCreatePointers: Oxy pnters;" );
490 for(
long i=1; i <= 8; i++ )
494 fprintf(
ioQQQ,
"\n" );
523 for(
long i=1; i <=
nLevel1; i++ )
545 else if(
TauLines[i].Emis().Aul() > 0. )
555 fprintf(
ioQQQ,
" level 1 line does not have valid gf or A\n" );
556 fprintf(
ioQQQ,
" This is ContCreatePointers\n" );
576 for (
int ipSpecies=0; ipSpecies <
nSpecies; ++ipSpecies)
579 em !=
dBaseTrans[ipSpecies].Emis().end(); ++em)
582 (*em).dampXvel() = (
realnum)(1./
583 dBaseStates[ipSpecies][em->Tran().ipHi()].lifetime()/em->Tran().EnergyWN()/
PI4);
584 (*em).damp() = -1000.0;
587 strncpy(chLab,(*(*em).Tran().Hi()).chLabel(),4);
590 static const double minAul = 1e-29;
591 if( (*em).Aul() > minAul )
593 (*em).Tran().ipCont() =
ipLineEnergy((*em).Tran().EnergyRyd(), chLab ,0);
594 (*em).ipFine() =
ipFineCont((*em).Tran().EnergyRyd() );
598 (*em).Tran().ipCont() = -1;
604 (*(*em).Tran().Lo()).g()));
611 (*diatom)->H2_ContPoint();
616 for(
long nelem=2; nelem <
LIMELM; nelem++ )
618 if(
dense.lgElmtOn[nelem])
621 sprintf( chLab,
"%2s%2ld",
elementnames.chElementSym[nelem], nelem+1-ipISO );
623 iso_sp[ipISO][nelem].fb[0].ipIsoLevNIonCon =
626 for(
long ipHi=1; ipHi <
iso_sp[ipISO][nelem].numLevels_max; ipHi++ )
632 for(
long ipLo=0; ipLo < ipHi; ipLo++ )
635 if(
iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().Aul() <=
iso_ctrl.SmallA )
642 if(
iso_sp[ipISO][nelem].trans(ipHi,ipLo).EnergyRyd() <
continuum.filbnd[0] )
645 iso_sp[ipISO][nelem].trans(ipHi,ipLo).ipCont() =
647 iso_sp[ipISO][nelem].fb[ipLo].ipIsoLevNIonCon);
648 iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().ipFine() =
659 for(
long nelem=ipISO; nelem <
LIMELM; nelem++ )
661 if(
dense.lgElmtOn[nelem])
664 for(
long ipHi=2; ipHi <
iso_ctrl.nLyman_malloc[ipISO]; ipHi++ )
672 iso_sp[ipISO][nelem].fb[ipLo].ipIsoLevNIonCon);
674 (*tr).Emis().ipFine() =
682 for(
long ipLo=0; ipLo<
iso_sp[ipISO][nelem].numLevels_max; ipLo++ )
702 for(
long nelem=ipISO; nelem <
LIMELM; nelem++ )
704 if(
dense.lgElmtOn[nelem])
706 for(
long ipHi=1; ipHi <
iso_sp[ipISO][nelem].numLevels_max; ipHi++ )
708 for(
long ipLo=0; ipLo < ipHi; ipLo++ )
727 for(
long i=0; i<
nUTA; ++i )
792 HFLines[i].Emis().damp() = 1e-20f;
831 for(
long nelem=ipISO; nelem<
LIMELM; ++nelem )
833 if(
dense.lgElmtOn[nelem] )
837 const int TwoS = (1+ipISO);
841 Aul = 8.226*pow((
double)(nelem+1.),6.);
849 const double As2nuFrom1S[29] = {51.02,1940.,1.82E+04,9.21E+04,3.30E+05,9.44E+05,2.31E+06,5.03E+06,1.00E+07,
850 1.86E+07,3.25E+07,5.42E+07,8.69E+07,1.34E+08,2.02E+08,2.96E+08,4.23E+08,5.93E+08,8.16E+08,
851 1.08E+09,1.43E+09,1.88E+09,2.43E+09,3.25E+09,3.95E+09,4.96E+09,6.52E+09,7.62E+09,9.94E+09};
852 Aul = As2nuFrom1S[nelem-1];
857 iso_sp[ipISO][nelem].trans(TwoS,0),
866 for(
long nelem=ipISO; nelem<
LIMELM; ++nelem )
868 if(
dense.lgElmtOn[nelem] )
876 const double As2nuFrom3S[29] = {4.09e-9,1.25E-06,5.53E-05,8.93E-04,8.05E-03,4.95E-02,2.33E-01,8.94E-01,2.95E+00,
877 8.59E+00,2.26E+01,5.49E+01,1.24E+02,2.64E+02,5.33E+02,1.03E+03,1.91E+03,3.41E+03,5.91E+03,
878 9.20E+03,1.50E+04,2.39E+04,3.72E+04,6.27E+04,8.57E+04,1.27E+05,2.04E+05,2.66E+05,4.17E+05};
881 As2nuFrom3S[nelem-1],
890 enum {DEBUG_LOC=
false};
895 0., 0.03738, 0.07506, 0.1124, 0.1498, 0.1875,
896 0.225, 0.263, 0.300, 0.3373, 0.375, 0.4127,
897 0.4500, 0.487, 0.525, 0.5625, 0.6002, 0.6376,
898 0.6749, 0.7126, 0.75};
906 for(
long i=0; i < nCRS; i++ )
908 fprintf(
ioQQQ,
"%.3e\t%.3e\n", ener[i] ,
914 for(
long i=0; i < limit; i++ )
917 fprintf(
ioQQQ,
"%.3e\t%.3e\t%.3e\n",
921 xnew += tnu.
As2nu[i];
923 fprintf(
ioQQQ,
" sum is %.3e\n", xnew );
929 enum {DEBUG_LOC=
false};
932 for(
long i=0; i<11; ++i )
935 (*TauDummy).WLAng() = (
realnum)(
PI * pow(10.,(
double)i));
937 fprintf(
ioQQQ,
"%.2f\t%s\n", (*TauDummy).WLAng() , chLsav );
946 fprintf(
ioQQQ,
" WL(Ang) E(RYD) IP gl gu gf A damp abs K\n" );
947 for(
long i=1; i <=
nLevel1; i++ )
950 long iWL_Ang = (long)
TauLines[i].WLAng();
951 if( iWL_Ang > 1000000 )
955 else if( iWL_Ang > 10000 )
960 fprintf(
ioQQQ,
" %10.10s%5ld%10.3e %4li%4ld%4ld%10.2e%10.2e%10.2e%10.2e\n",
969 for (
int ipSpecies=0; ipSpecies <
nSpecies; ++ipSpecies)
972 em !=
dBaseTrans[ipSpecies].Emis().end(); ++em)
976 long iWL_Ang = (long)(*em).Tran().WLAng();
978 if( iWL_Ang > 1000000 )
982 else if( iWL_Ang > 10000 )
986 fprintf(
ioQQQ,
" %10.10s%5ld%10.3e %4li%4ld%4ld%10.2e%10.2e%10.2e%10.2e\n",
987 chLab, iWL_Ang,
RYDLAM/(*em).Tran().WLAng(),
988 (*em).Tran().ipCont(), (
long)((*(*em).Tran().Lo()).g()),
989 (
long)((*(*em).Tran().Hi()).g()),(*em).gf(),
990 (*em).Aul(),(*em).dampXvel(),
999 long iWL_Ang = (long)
TauLine2[i].WLAng();
1001 if( iWL_Ang > 1000000 )
1005 else if( iWL_Ang > 10000 )
1009 fprintf(
ioQQQ,
" %10.10s%5ld%10.3e %4li%4ld%4ld%10.2e%10.2e%10.2e%10.2e\n",
1020 long iWL_Ang = (long)
HFLines[i].WLAng();
1022 if( iWL_Ang > 1000000 )
1026 else if( iWL_Ang > 10000 )
1030 fprintf(
ioQQQ,
" %10.10s%5ld%10.3e %4li%4ld%4ld%10.2e%10.2e%10.2e%10.2e\n",
1035 HFLines[i].Emis().opacity() );
1042 for(
long i=1; i <=
nLevel1; i++ )
1044 if(
TauLines[i].EnergyWN() < 10000. )
1051 for (
int ipSpecies=0; ipSpecies <
nSpecies; ++ipSpecies)
1054 em !=
dBaseTrans[ipSpecies].Emis().end(); ++em)
1056 if((*em).Tran().EnergyWN() < 10000. )
1058 (*em).opacity() = 0.;
1138 double thresh=-DBL_MAX;
1154 for( ion=0; ion < nelem; ion++ )
1164 long int ipISO = nelem-ion;
1171 imax =
Heavy.nsShells[nelem][ion];
1174 for( nshell=0; nshell < imax; nshell++ )
1184 opac.ipElement[nelem][ion][nshell][0] = 2;
1185 opac.ipElement[nelem][ion][nshell][1] = 1;
1193 opac.ipElement[nelem][ion][nshell][0] =
1203 opac.ipElement[nelem][ion][nshell][1] =
1204 LimitSh(ion+1, nshell+1,nelem+1);
1205 ASSERT(
opac.ipElement[nelem][ion][nshell][1] > 0);
1209 ASSERT( imax > 0 && imax <= 7 );
1213 opac.ipElement[nelem][ion][imax-1][0] =
1217 Heavy.ipHeavy[nelem][ion] =
opac.ipElement[nelem][ion][imax-1][0];
1222 Heavy.Valence_IP_Ryd[nelem][ion] = thresh;
1224 Heavy.xLyaHeavy[nelem][ion] = 0.;
1229 Heavy.ipLyHeavy[nelem][ion] =
1232 Heavy.ipBalHeavy[nelem][ion] =
1239 Heavy.ipLyHeavy[nelem][ion] = -1;
1240 Heavy.ipBalHeavy[nelem][ion] = -1;
1246 Heavy.nsShells[nelem][nelem] = 1;
1255 ASSERT(
opac.ipElement[nelem][nelem][0][0] > 0 );
1260 Heavy.ipHeavy[nelem][nelem] =
opac.ipElement[nelem][nelem][0][0];
1265 for( ion=0; ion < (nelem+1); ion++ )
1267 fprintf(
ioQQQ,
"Ion:%3ld%3ld %2.2s%2.2s total shells:%3ld\n",
1269 ,
Heavy.nsShells[nelem][ion] );
1270 for( ns=0; ns <
Heavy.nsShells[nelem][ion]; ns++ )
1272 fprintf(
ioQQQ,
" shell%3ld %2.2s range eV%10.2e-%8.2e\n",
1273 ns+1,
Heavy.chShell[ns],
rfield.anu[
opac.ipElement[nelem][ion][ns][0]-1]*