/* This file is part of Cloudy and is copyright (C)1978-2006 by Gary J. Ferland
 * For conditions of distribution and use see copyright notice in license.txt */
/* routine called at end of iteration when advection is turned on */
void DynaEndIter(void);

/* DynaStartZone called at start of iteration when advection is turned on */
void DynaStartZone(void);

/* DynaEndZone called at end of iteration when advection is turned on */
void DynaEndZone(void);

/* DynaIonize, called from ionize to evaluate advective terms for current conditions */
void DynaIonize(void);

/* DynaPresChngFactor, called from PressureChange to evaluate new density needed for
 * current conditions and wind solution, returns ratio of new to old density,
 * called when wind velocity is negative for D-critical solutions  */
double DynaPresChngFactor(void);

/*DynaZero zero some dynamics variables, called from zero.c */
void DynaZero( void );

/*DynaCreateArrays allocate some space needed to save the dynamics structure variables, called from atmdat_readin */
void DynaCreateArrays( void );

/* ParseDynaWind parse the wind command, called from ParseCommands */
void ParseDynaWind( char *chCard );

/* DynaPrtZone - called to print zone results */
void DynaPrtZone( void );

/*DynaPunch punch info related to advection */
void DynaPunch(FILE* ipPnunit , char chJob );

float DynaFlux(double depth);

/* all of these are initialized in zero */
EXTERN struct t_advec 
{
	/* is advection turned on ?, set to FALSE in zero */
	int lgAdvection;

	/* advective cooling minus heating (erg cm^-3 s^-1) */
	double Cool, Heat, dCooldT, dHeatdT;

	/* largest fraction of cooling and heating */
	double CoolMax, HeatMax;

	/* the advection rate (s^-1) */
	double Rate;

	/* the advective recombination rate (s^-1)*/
	double **Source /*[LIMELM][LIMELM+2]*/;

	double Molec[N_H_MOLEC];

	/* factor to turn off advection for H-like or He-like iso seq, no advection h-like, he-like*/
	int lgISO[NISO];
	/* factor to turn off advection for He-like iso seq, no advection h-like*/
	int lgHeLIKE;
	/* factor to turn off advection for rest of ions, treated in bidiag, no advection h-like*/
	int lgMETALS;
	/* factor to turn off advective cooling */
	int lgCoolHeat;

	/* the initial value of the advection length, reset with set advection length */
	double AdvecLengthInit;

 	/* the center of the particle flux law */
	double FluxCenter;

	/* flag set by the "set dynamics pressure mode" command */
	char chPresMode[20];

	/* the shock depth in cm set with "set dynamics shock depth" command */
	double ShockDepth;

	/* the isothermal Mach number at which to insert an antishock
	   set with "set dynamics antishock mach" command */
	double ShockMach;

 	/* the scale of the particle flux law */
	double FluxScale;

	/* whether we also need to scale by the face density */
	int lgFluxDScale;

 	/* the power law index of the particle flux law */
	double FluxIndex;

	/* the proposed thickness for the next zone when advection is included */
	double dRad;

	/* the depth of the last iteration */
	double oldFullDepth;

	/* convergence_error and discretization_error give estimates of convergence: 
	 :: convergence_error -- change between the last iterations;
	 :: discretization_error -- error in the upstream interpolation.
	 When (and if) discretization_error >> convergence_error, the interpolation length should be decreased.

	 They should both be based on the same norm of the models, but what norm may be
	 experimented with -- at present, it's H+/Htot just weighted by cell number, 
	 which makes the estimates sensitive to the structure of the primary ionization front.
	 */
	/* the error from comparing this iteration with the previous one */
	double convergence_error ;

	/* the allowed rel error, by default 0.1 */
	double convergence_tolerance;

	/* the error to be expected based on the coarseness of current advection length */
	double discretization_error ;

	/* two ways of scaling the error estimate for convergence */
	double error_scale1 , error_scale2 ;

	/* flag set true if set dynamics flow type was set - this means
	 * to use the specified option, and not to derive one */
	int lgSetPresMode;

	/* deriv of density needed for expansion cooling */
	float dDensityDT;
	/* first iteration where dynamics is included - default is 2 - the second true iteration,
	 * set with SET DYNAMICS FIRST command */
	long iteration_first_dynamics;

}	dynamics;
