Merge branch '22apr' into 'master'

Protection against division by zero cooling.

See merge request !35

source/age_check.cpp

@@ -39,14 +39,17 @@ void AgeCheck(void)
 
 	for( i=0; i < limit; i++ )
 	{
-		timesc.time_therm_long = 
-			MAX2( timesc.time_therm_long ,
-			struc.DenParticles[i]*BOLTZMANN*1.5*struc.testr[i]/struc.coolstr[i]);
-		timesc.time_therm_short = 
-			MIN2( timesc.time_therm_short ,
-			struc.DenParticles[i]*BOLTZMANN*1.5*struc.testr[i]/struc.coolstr[i]);
-		/*>>chng 99 feb 01, had div by heating, changed to cooling so constant
-		 * temperature models are more realistic */
+		if( struc.coolstr[i]>SMALLFLOAT )
+		{
+			timesc.time_therm_long =
+				MAX2( timesc.time_therm_long ,
+				struc.DenParticles[i]*BOLTZMANN*1.5*struc.testr[i]/struc.coolstr[i]);
+			timesc.time_therm_short =
+				MIN2( timesc.time_therm_short ,
+				struc.DenParticles[i]*BOLTZMANN*1.5*struc.testr[i]/struc.coolstr[i]);
+			/*>>chng 99 feb 01, had div by heating, changed to cooling so constant
+			 * temperature models are more realistic */
+		}
 	}
 
 	tlong = MAX2(tlong,timesc.time_therm_long);

source/cool_save.cpp

@@ -67,7 +67,10 @@ void CoolSave(FILE * io, const char chJob[])
 		if( strncmp( thermal.chClntLab[i], "adve", 4 ) == 0 )
 			continue;
 
-		csav[ip] = (realnum)( safe_div( MAX2(thermal.cooling[i],thermal.heatnt[i]), cool_total, 0. ));
+		if( cool_total>SMALLFLOAT )
+			csav[ip] = (realnum)( safe_div( MAX2(thermal.cooling[i],thermal.heatnt[i]), cool_total, 0. ));
+		else
+			csav[ip] = 0.;
 		/* save sign to remember if heating or cooling line */
 		if( thermal.heatnt[i] == 0. )
 		{

source/dense.cpp

@@ -179,6 +179,9 @@ bool lgElemsConserved (void)
 				sum_monatomic + sum_molecules,
 				sum_gas_phase,
 				(sum_monatomic + sum_molecules - sum_gas_phase)/sum_gas_phase );
+			fprintf(ioQQQ," Consider increasing value of conv.GasPhaseAbundErrorAllowed with SET DENSITY TOLERANCE as a workaround. "
+					"Its current value is %.2e\n" , conv.GasPhaseAbundErrorAllowed );
+			fprintf(ioQQQ," Temperature %.2e H density %.2e zone %ld \n", phycon.te, dense.gas_phase[ipHYDROGEN], nzone);
 			lgOK=false;
 		}
 		

source/prt_comment.cpp

@@ -792,35 +792,35 @@ void PrtComment(void)
 	if( !rfield.lgMMok )
 	{
 		sprintf( chLine, 
-			" C-Continuum not defined in extreme infrared - Compton scat, grain heating, not treated properly?" );
-		caunin(chLine);
+			"  !Continuum not defined in extreme infrared - Compton scat, grain heating, not treated properly?" );
+		bangin(chLine);
 	}
 
 	if( !rfield.lgHPhtOK )
 	{
 		sprintf( chLine, 
-			" C-Continuum not defined at photon energies which ionize excited states of H, important for H- and ff heating." );
-		caunin(chLine);
+			"  !Continuum not defined at photon energies which ionize excited states of H, important for H- and ff heating." );
+		bangin(chLine);
 	}
 
 	if( !rfield.lgXRayOK )
 	{
 		sprintf( chLine, 
-			" C-Continuum not defined at X-Ray energies - Compton scattering and Auger ionization wrong?" );
-		caunin(chLine);
+			"  !Continuum not defined at X-Ray energies - Compton scattering and Auger ionization wrong?" );
+		bangin(chLine);
 	}
 
 	if( !rfield.lgGamrOK )
 	{
 		sprintf( chLine, 
-			" C-Continuum not defined at gamma-ray energies - pair production and Compton scattering OK?" );
-		caunin(chLine);
+			"  !Continuum not defined at gamma-ray energies - pair production and Compton scattering OK?" );
+		bangin(chLine);
 	}
 
 	if( continuum.lgCon0 )
 	{
-		sprintf( chLine, " C-Continuum zero at some energies." );
-		caunin(chLine);
+		sprintf( chLine, "  !Continuum zero at some energies." );
+		bangin(chLine);
 	}
 
 	if( continuum.lgCoStarInterpolationCaution )

source/timesc.cpp

@@ -41,12 +41,16 @@ void t_timesc::calc_therm_timesc( long int izone )
 	{
 		/* >>chng 99 feb 01, had div by heating, changed to cooling so constant
 		 * temperature models are more realistic */
-		double dt = 1.5 * BOLTZMANN * struc.DenParticles[i] * struc.testr[i] / struc.coolstr[i];
-		time_therm_long = MAX2( time_therm_long , dt );
-		time_therm_short= MIN2( time_therm_short, dt );
-		//      printf("dt = %g\t long = %g\t short = %g\n", dt, time_therm_long, time_therm_short);
-	}
+		if( struc.coolstr[i] > SMALLFLOAT )
+		{
+			double dt = 1.5 * BOLTZMANN * struc.DenParticles[i] * struc.testr[i] / struc.coolstr[i];
+			//double dt = 1.5 * BOLTZMANN * struc.DenParticles[i] * safe_div( (double)struc.testr[i] , struc.coolstr[i] );
+			time_therm_long = MAX2( time_therm_long , dt );
+			time_therm_short= MIN2( time_therm_short, dt );
+			//      printf("dt = %g\t long = %g\t short = %g\n", dt, time_therm_long, time_therm_short);}
+		}
 	//      printf( "*** long = %g\t short = %g\n", time_therm_long, time_therm_short );
+	}
 
 	return;
 }