/*!\file:  IsInputConverged.cpp
 * \brief: figure out convergence at the input level. 
 * We are given a list of inputs, new and old, and a criterion, and using the inputs, we 
 * return the value of the criterion test, which will be used at the solution level 
 * to determine convergence.
 */ 

#ifdef HAVE_CONFIG_H
	#include <config.h>
#else
#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
#endif

#include "../../classes/objects/objects.h"

void IsInputConverged(IssmDouble* peps, Input** new_inputs,Input** old_inputs,int num_inputs,int criterion_enum){

	int i,j;

	/*output: */
	IssmDouble eps;
	
	/*intermediary: */
	IssmDouble *newvalues     = NULL;
	IssmDouble *oldvalues     = NULL;
	int     num_values;
	IssmDouble  ndu        = 0;
	IssmDouble  nu         = 0;

	if(criterion_enum==RelativeEnum){

		/*conpute ndu/du (where u could be velocity, pressure, temperature, etc ...): */
		for(i=0;i<num_inputs;i++){

			/*in debugging mode, check that the inputs are of the same type*/
			_assert_(new_inputs[i]->ObjectEnum()==old_inputs[i]->ObjectEnum());

			/*Get pointers*/
			new_inputs[i]->GetValuesPtr(&newvalues,&num_values);
			old_inputs[i]->GetValuesPtr(&oldvalues,&num_values);
			for(j=0;j<num_values;j++){
				ndu+=pow(newvalues[j]-oldvalues[j],2);
				nu+=pow(oldvalues[j],2);
			}
		}

		/*take square root: */
		ndu=sqrt(ndu);
		nu=sqrt(nu);

		/*now, compute eps: */
		if(reCast<bool>(nu))eps=ndu/nu;
		else eps=0;
	}
	else _error2_("convergence criterion " << EnumToStringx(criterion_enum) << " not supported yet!");

	/*Assign output pointers:*/
	*peps=eps;
}
