From 1d21f0b65dba365bb5d4d1788992247bd56fc786 Mon Sep 17 00:00:00 2001 From: Baptiste Mouginot <mouginot.baptiste@gmail.com> Date: Thu, 2 Jul 2015 10:37:00 +0000 Subject: [PATCH] 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 | 11 +--- source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx | 4 +- .../BaM_Dev/include/EquivalenceModel.hxx | 4 ++ .../branches/BaM_Dev/src/EquivalenceModel.cxx | 61 ++++++++++++++++--- 4 files changed, 63 insertions(+), 17 deletions(-) diff --git a/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx b/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx index 48f96a73d..8d7667be8 100644 --- a/source/branches/BaM_Dev/Model/XS/XSM_MLP.cxx +++ b/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(); diff --git a/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx b/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx index 22ab951bb..10001e310 100644 --- a/source/branches/BaM_Dev/Model/XS/XSM_MLP.hxx +++ b/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(); diff --git a/source/branches/BaM_Dev/include/EquivalenceModel.hxx b/source/branches/BaM_Dev/include/EquivalenceModel.hxx index 2683cc865..deaf6be75 100644 --- a/source/branches/BaM_Dev/include/EquivalenceModel.hxx +++ b/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 diff --git a/source/branches/BaM_Dev/src/EquivalenceModel.cxx b/source/branches/BaM_Dev/src/EquivalenceModel.cxx index c05778266..6769842e8 100644 --- a/source/branches/BaM_Dev/src/EquivalenceModel.cxx +++ b/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; -} \ No newline at end of file + for(int i=0 ; i < (int)lambda.size() ;i++ ) + lambda[i]= -1; +} + + + +//________________________________________________________________________ +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; + +} + -- GitLab