Add isInDomain to EQ model

git-svn-id: svn+ssh://svn.in2p3.fr/class@685 0e7d625b-0364-4367-a6be-d5be4a48d228

source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx

@@ -39,10 +39,7 @@ XSM_MLP::XSM_MLP(string TMVA_Weight_Directory,string InformationFile, bool IsTim
 	fIsStepTime = IsTimeStep;
 	fTMVAWeightFolder = TMVA_Weight_Directory;
 	
-	if(InformationFile=="")
-		fInformationFile = TMVA_Weight_Directory + "/Data_Base_Info.nfo";
-	else
-		fInformationFile = fTMVAWeightFolder+InformationFile;
+	fInformationFile = fTMVAWeightFolder+InformationFile;
 	
 	GetMLPWeightFiles();
 	
@@ -61,10 +58,8 @@ XSM_MLP::XSM_MLP(CLASSLogger* Log,string TMVA_Weight_Directory,string Informatio
 	
 	fIsStepTime = IsTimeStep;
 	fTMVAWeightFolder = TMVA_Weight_Directory;
-	if( InformationFile == "" )
-		fInformationFile = TMVA_Weight_Directory + "/Data_Base_Info.nfo";
-	else
-		fInformationFile = fTMVAWeightFolder+InformationFile;
+
+	fInformationFile = fTMVAWeightFolder+InformationFile;
 	
 	GetMLPWeightFiles();
 	

source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx

@@ -60,7 +60,7 @@ class XSM_MLP : public XSModel
 	 \param IsTimeStep : if true , one TMVA weihgt per step time is requiered otherwise it assumes time is part of the MLP inputs
 	 
 	 */
-	XSM_MLP(string TMVA_Weight_Directory,string InformationFile="",bool IsTimeStep=false);
+	XSM_MLP(string TMVA_Weight_Directory,string InformationFile="/Data_Base_Info.nfo",bool IsTimeStep=false);
 	//}
 	
 	//{
@@ -72,7 +72,7 @@ class XSM_MLP : public XSModel
 	 \param IsTimeStep : if true , one TMVA weihgt per step time is requiered otherwise it assumes time is part of the MLP inputs
 	 
 	 */
-	XSM_MLP(CLASSLogger* Log,string TMVA_Weight_Directory,string InformationFile="",bool IsTimeStep=false);
+	XSM_MLP(CLASSLogger* Log,string TMVA_Weight_Directory,string InformationFile="/Data_Base_Info.nfo",bool IsTimeStep=false);
 	//}
 	
 	~XSM_MLP();

source/branches/BaM_Dev/include/EquivalenceModel.hxx

@@ -110,6 +110,10 @@ class EquivalenceModel : public CLASSObject
 	void ReadType(const string &line);
 	//@}
 	
+	
+	bool isIVInDomain(IsotopicVector IV);
+
+	
 	protected :
 	
 	IsotopicVector fFertileList;	//!< contain the list of zai, needed as fertile, taken in a stock before fabrication

source/branches/BaM_Dev/src/EquivalenceModel.cxx

@@ -8,7 +8,7 @@ EquivalenceModel::EquivalenceModel():CLASSObject()
 	fMaxInterration = 100; // Max iterration in build fueld algorythum
 	fFirstGuessFissilContent = 0.02;
 	freaded = false;
-
+	
 }
 
 EquivalenceModel::EquivalenceModel(CLASSLogger* log):CLASSObject(log)
@@ -259,7 +259,7 @@ vector<double> EquivalenceModel::BuildFuel(double BurnUp, double HMMass, vector<
 		if( LAMBDA_NEEDED == -1 )	// Check if previous lambda was well calculated
 		{
 			SetLambdaToErrorCode(lambda);
-			WARNING << "Not enought fissile material to build fuel" << endl;		
+			WARNING << "Not enought fissile material to build fuel" << endl;
 			return lambda;
 		}
 		
@@ -290,11 +290,11 @@ vector<double> EquivalenceModel::BuildFuel(double BurnUp, double HMMass, vector<
 			WARNING<<"GetFissileMolarFraction return negative or greater than one value";
 			return lambda;
 		}
-
+		
 		double MeanMolarPu = Fissile.GetMeanMolarMass();
 		double MeanMolarDepletedU = Fertile.GetMeanMolarMass();
 		
-		double MeanMolar = MeanMolarPu * MolarPuContent + (1-MolarPuContent)  *MeanMolarDepletedU;
+		double MeanMolar = MeanMolarPu * MolarPuContent + (1-MolarPuContent)  * MeanMolarDepletedU;
 		
 		
 		WeightPuContent = MolarPuContent * MeanMolarPu / MeanMolar;
@@ -307,9 +307,11 @@ vector<double> EquivalenceModel::BuildFuel(double BurnUp, double HMMass, vector<
 		
 	}while(  fabs( PuMassNeeded - AvailablePuMass )/HMMass > fRelativMassPrecision );
 	
+	(*this).isIVInDomain(fissil);
 	
 	DBGV( "Weight percent fissil : " << PuMassNeeded/HMMass );
 	DBGV( "Lambda vector : " );
+	
 	for(int i = 0; i < (int)FissilArray.size() + (int)FertilArray.size(); i++ )
 		DBGV(lambda[i]);
 	
@@ -335,9 +337,54 @@ void EquivalenceModel::SetLambda(vector<double>& lambda ,int FirstStockID, int L
 	
 	lambda[FirstStockID + IntegerPart] = DecimalPart;
 }
+
 //________________________________________________________________________
 void EquivalenceModel::SetLambdaToErrorCode(vector<double>& lambda)
 {
-		for(int i=0 ; i < (int)lambda.size() ;i++ )
-			lambda[i]= -1;
-}
+	for(int i=0 ; i < (int)lambda.size() ;i++ )
+		lambda[i]= -1;
+}
\ No newline at end of file
+
+
+
+//________________________________________________________________________
+bool EquivalenceModel::isIVInDomain(IsotopicVector IV)
+{
+	DBGL
+	bool IsInDomain=true;
+	
+	if(fZAILimits.empty())
+	{
+	 WARNING << "Fresh Fuel variation domain is not set" << endl;
+	 WARNING << "CLASS has no clue if the computed evolution for this fresh fuel is correct" << endl;
+	 WARNING << "Proceed finger crossed !!" << endl;
+	 return true;
+	}
+	
+	else
+	{
+		IsotopicVector IVNorm = IV /IV.GetSumOfAll();
+		for (map< ZAI,pair<double,double> >::iterator Domain_it=fZAILimits.begin(); Domain_it!=fZAILimits.end(); Domain_it++)
+		{
+			double ThatZAIProp = IVNorm.GetIsotopicQuantity()[Domain_it->first]	;
+			double ThatZAIMin  = Domain_it->second.first;
+			double ThatZAIMax  = Domain_it->second.second;
+			if( (ThatZAIProp > ThatZAIMax) || (ThatZAIProp <  ThatZAIMin) )
+			{
+				IsInDomain = false;
+				
+				WARNING << "Fresh fuel out of model range" << endl;
+				WARNING << "\t AT LEAST this ZAI is accused to be outrange :" << endl;
+				WARNING << "\t\t" << Domain_it->first.Z() << " "<<Domain_it->first.A() << " " << Domain_it->first.I() << endl;
+				WARNING << "\t\t min=" << ThatZAIMin <<" value=" << ThatZAIProp << " max=" << ThatZAIMax << endl;
+				WARNING << "\t IV accused :" << endl << endl;
+				WARNING << IVNorm.sPrint() << endl;
+				break;
+			}
+		}
+	}
+	DBGL
+	return IsInDomain;
+	
+}
+