/* This file is part of Cloudy and is copyright (C) 1978-2004 by Gary J. Ferland.
 * For conditions of distribution and use, see copyright notice in license.txt */
/*CoolSum  total cooling from all entries into cooling stack */
/*CoolZero set cooling and heating stack to zero */
/*CoolAdd add coolants to the cooling stack, called in evaluation of cooling function */
#include "cddefines.h"
#include "taulines.h"
#include "lines_service.h"
#include "thermal.h"
#include "cooling.h"

/*CoolAdd add coolants to the cooling stack, called in evaluation of cooling function */
void CoolAdd(char *chLabel, 
  float lambda, 
  double cool)
{

#	ifdef DEBUG_FUN
	fputs( "<+>CoolAdd()\n", debug_fp );
#	endif

	/* this flag indicates (TRUE) that we are between when cooling was set to
	 * zero with call to CoolZero, and when final sum was used.  Any call
	 * after final summation (FALSE) will be ignored and so is fatal error */
	ASSERT( thermal.lgCoolEvalOK );

	/* this can be done with an assert since these results cannot possibly
	 * depend on user input */
	ASSERT( thermal.ncltot < NCOLNT );

	/* copy coolant label into stack */
	strcpy( thermal.chClntLab[thermal.ncltot], chLabel);

	/* now the wavelength */
	thermal.collam[thermal.ncltot] = lambda;

	/* normal line cooling */
	thermal.cooling[thermal.ncltot] = MAX2(0.,cool);

	/* possible line heating - not supposed to be done this way!
	 * this is intrinsic positive number, to be added to heating */
	thermal.heatnt[thermal.ncltot] = MAX2(0.,-cool);

	/* now increment counter, this is the number of coolants entered */
	thermal.ncltot += 1;

#	ifdef DEBUG_FUN
	fputs( " <->CoolAdd()\n", debug_fp );
#	endif
	return;
}

/*CoolZero set cooling and heating stack to zero */
void CoolZero(void)
{

#	ifdef DEBUG_FUN
	fputs( "<+>CoolZero()\n", debug_fp );
#	endif

	thermal.ncltot = 0;
	thermal.dCooldT = 0.;

	/* >>chng 03 nov 29, from explicit loop to memset to save time */
	memset(thermal.cooling , 0 , NCOLNT*sizeof(thermal.cooling[0] ) );
	memset(thermal.heatnt , 0 , NCOLNT*sizeof(thermal.heatnt[0] ) );
#	if 0
	for( i=0; i < NCOLNT; i++ )
	{
		thermal.cooling[i] = 0.;
		thermal.heatnt[i] = 0.;
	}
#	endif

	/* this flag indicates that it is ok to add coolants to cooling
	 * stack since between first zero, and final sum - CoolAdd checks
	 * that this is TRUE */
	thermal.lgCoolEvalOK = TRUE;

#	ifdef DEBUG_FUN
	fputs( " <->CoolZero()\n", debug_fp );
#	endif
	return;
}

/*CoolSum  total cooling from all entries into cooling stack */
void CoolSum(double *total)
{
	long int i;

#	ifdef DEBUG_FUN
	fputs( "<+>CoolSum()\n", debug_fp );
#	endif

	/* routine to add together all line heating and cooling */

	*total = 0.;
	thermal.coolheat = 0.;
	/* this is sum of agents that should be coolants
	 * coolheat will be coolants that came out as heat */
	for( i=0; i < thermal.ncltot; i++ )
	{
		*total += thermal.cooling[i];
		thermal.coolheat += thermal.heatnt[i];
	}
	thermal.heating[0][12] = thermal.coolheat;

	/* make comment if negative cooling ever significant */
	if( thermal.htot > 0. )
	{
		if( thermal.coolheat/thermal.htot > 0.01 )
		{
			/* CoolHeatMax was set to zero at start of calc, we want very biggest */
			for( i=0; i < thermal.ncltot; i++ )
			{
				if( thermal.heatnt[i]/thermal.htot > thermal.CoolHeatMax )
				{
					thermal.CoolHeatMax = (float)(thermal.heatnt[i]/thermal.htot);
					thermal.wlCoolHeatMax = thermal.collam[i];
					strcpy( thermal.chCoolHeatMax, thermal.chClntLab[i] );
				}
			}
		}
	}

	/* this sum of lines that were heat sources - this
	 * part was not counted as heating in call to cooling add routine
	 * since atom_level2 and atom_level3 cooling routines separate this out 
	 * into ->cool and ->heat - this does
	 * NOT double count line heating */

	thermal.heatl = 0.;
	for( i=0; i < nWindLine; i++ )
	{
		if( TauLine2[i].IonStg < TauLine2[i].nelem+1-NISO )
			thermal.heatl += TauLine2[i].heat;
	}

	for( i=1; i <= nLevel1; i++ )
	{
		thermal.heatl += TauLines[i].heat;
	}
#	if 0
	/* >>chng 03 feb 25, this should not be counted here since
	 * total cooling (less heating) was added into cooling stack */
	for( i=0; i < nCORotate; i++ )
	{
		thermal.heatl += C12O16Rotate[i].heat;
		thermal.heatl += C13O16Rotate[i].heat;
	}
#	endif

	/* line heating added in following, only here */
	thermal.heating[0][22] = thermal.heatl;
	{
		/*@-redef@*/
		enum {DEBUG_LOC=FALSE};
		/*@+redef@*/
		if( DEBUG_LOC && thermal.heatl/thermal.ctot > 0.1 )
		{
			double thresh = 0.1;
			fprintf(ioQQQ," all heating lines > %.4f of heatl printed next \n",
				thresh );
			for( i=0; i < nWindLine; i++ )
			{
				if( TauLine2[i].IonStg < TauLine2[i].nelem+1-NISO )
				{
					if( TauLine2[i].heat/thermal.heatl > thresh )
						DumpLine( &TauLine2[i] );
				}
			}

			for( i=1; i <= nLevel1; i++ )
			{
				if( TauLines[i].heat/thermal.heatl > thresh )
					DumpLine( &TauLines[i] );
			}
		}
	}

	/*begin sanity check */
	if( *total <= 0. )
	{
		fprintf( ioQQQ, " CoolSum finds cooling <= 0%10.2e\n", 
		  *total );
	}
	if( thermal.heatl/thermal.ctot < -1e-15 )
	{
		fprintf( ioQQQ, " CoolSum finds negative heating %10.2e %10.2e\n", 
		  thermal.heatl, thermal.ctot );
		/* write(qq,'('' dump went to fort.66'')') 
		 * write(66,'(1p,12e10.1)') (heatnt(i),i=1,ncltot) */
	}
	/*end sanity check */

	thermal.lgCoolEvalOK = FALSE;

#	ifdef DEBUG_FUN
	fputs( " <->CoolSum()\n", debug_fp );
#	endif
	return;
}

