diff --git a/DATA_BASES/PWR/MOX/EQModel/MLP_Kinf/MLP/PWR_MOX_30Wg.nfo b/DATA_BASES/PWR/MOX/EQModel/MLP_Kinf/MLP/PWR_MOX_30Wg.nfo index 1d09819ee4844d0c90de9e68bcc2f363cb905ce9..a4428bb946e90e940065935c3a40587ceda50812 100644 --- a/DATA_BASES/PWR/MOX/EQModel/MLP_Kinf/MLP/PWR_MOX_30Wg.nfo +++ b/DATA_BASES/PWR/MOX/EQModel/MLP_Kinf/MLP/PWR_MOX_30Wg.nfo @@ -8,10 +8,7 @@ Specific Power (W/gHM) : k_specpower 34.24 Maximal burnup (GWd/tHM) : -k_maxburnup 500 - -Maximal fissile content (molar proportion) : -k_maxfiscontent 0.25 +k_maxburnup 75 Z A I Name (input MLP) : k_zainame 92 235 0 U5 @@ -36,9 +33,11 @@ Fertile Liste (Z A I Default Proportion) : k_list Fertile 92 235 0 0.0025 k_list Fertile 92 238 0 0.9975 -Starting fissile content in fuel for equivalence model calculation : -k_firstguesscontent Fissile 0.04 -k_firstguesscontent Fertile 0.96 +Minimum Fraction of each list in the fuel (%mass) : +k_MassFractionMin Fissile 0.04 + +Maximun Fraction of each list in the fuel (%mass) : +k_MassFractionMax Fertile 0.16 Date : 2014 Author : BaL diff --git a/source/Model/Equivalence/EQM_MLP_Kinf.cxx b/source/Model/Equivalence/EQM_MLP_Kinf.cxx index d2154de28851dd79f1572a4165ede0e789430f95..46dfedb7e5950c039b81cdd5372fbd62aa6dc0e8 100644 --- a/source/Model/Equivalence/EQM_MLP_Kinf.cxx +++ b/source/Model/Equivalence/EQM_MLP_Kinf.cxx @@ -34,16 +34,16 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(string WeightPathAlpha0, string WeightPathAlpha1, str fTMVAWeightPath.push_back(WeightPathAlpha2); fMaximalBU = 0; - fMaximalContent = 0; fNumberOfBatch = NumOfBatch; fKThreshold = CriticalityThreshold ; - SetBurnUpPrecision(0.005);//1 % of the targeted burnup - fInformationFile = InformationFile; LoadKeyword(); ReadNFO();//Getting information from fInformationFile + SetBurnUpPrecision(0.008);//0.8 % of the targeted burnup + SetPCMPrecision(10); + INFO << "__An equivalence model has been define__" << endl; INFO << "\tThis model is based on the prediction of kinf" << endl; INFO << "\tThe TMVA weight files are :" << endl; @@ -52,7 +52,6 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(string WeightPathAlpha0, string WeightPathAlpha1, str INFO << "\t" << fTMVAWeightPath[2] << endl; INFO << "\tThe Information file is :" << endl; INFO << "\t" << fInformationFile << endl; - INFO << "Maximal fissile content (molar proportion) : " << fMaximalContent << endl; INFO << "Maximal burnup (GWd/tHM) : " << fMaximalBU << endl; EquivalenceModel::PrintInfo(); @@ -66,13 +65,7 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(string WeightPathAlpha0, string WeightPathAlpha1, str if(fStreamList[(* it_s_IV).first].GetIsotopicQuantity().empty()){EmptyList=true;} } - for( it_s_D = fFirstGuessContent .begin(); it_s_D != fFirstGuessContent.end(); it_s_D++) - { - if(fFirstGuessContent[(* it_s_D).first] ==0){FirstContentNULL=true;} - } - - - if(fMapOfTMVAVariableNames.empty() || EmptyList==true || FirstContentNULL==true || fMaximalBU == 0 || fMaximalContent == 0 ) + if(fMapOfTMVAVariableNames.empty() || EmptyList==true || fMaximalBU == 0 |) { ERROR<<"Missing information file in : "<<fInformationFile<<endl; exit(1); @@ -88,16 +81,16 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(CLASSLogger* log, string WeightPathAlpha0, string Wei fTMVAWeightPath.push_back(WeightPathAlpha2); fMaximalBU = 0; - fMaximalContent = 0; fNumberOfBatch = NumOfBatch; fKThreshold = CriticalityThreshold ; - - SetBurnUpPrecision(0.005);//1 % of the targeted burnup fInformationFile = InformationFile; LoadKeyword(); ReadNFO();//Getting information from fInformationFile + SetBurnUpPrecision(0.008);//0.8 % of the targeted burnup + SetPCMPrecision(10); + INFO << "__An equivalence model has been define__" << endl; INFO << "\tThis model is based on the prediction of kinf" << endl; INFO << "\tThe TMVA weight files are :" << endl; @@ -106,7 +99,6 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(CLASSLogger* log, string WeightPathAlpha0, string Wei INFO << "\t" << fTMVAWeightPath[2] << endl; INFO << "\tThe Information file is :" << endl; INFO << "\t" << fInformationFile << endl; - INFO << "Maximal fissile content (molar proportion) : " << fMaximalContent << endl; INFO << "Maximal burnup (GWd/tHM) : " << fMaximalBU << endl; EquivalenceModel::PrintInfo(); @@ -120,13 +112,7 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(CLASSLogger* log, string WeightPathAlpha0, string Wei if(fStreamList[(* it_s_IV).first].GetIsotopicQuantity().empty()){EmptyList=true;} } - for( it_s_D = fFirstGuessContent .begin(); it_s_D != fFirstGuessContent.end(); it_s_D++) - { - if(fFirstGuessContent[(* it_s_D).first] ==0){FirstContentNULL=true;} - } - - - if(fMapOfTMVAVariableNames.empty() || EmptyList==true || FirstContentNULL==true || fMaximalBU == 0 || fMaximalContent == 0 ) + if(fMapOfTMVAVariableNames.empty() || EmptyList==true || fMaximalBU == 0) { ERROR<<"Missing information file in : "<<fInformationFile<<endl; exit(1); @@ -145,11 +131,10 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(string TMVAWeightPath, int NumOfBatch, string Inform InformationFile = StringLine::ReplaceAll(TMVAWeightPath,".xml",".nfo"); fMaximalBU = 0; - fMaximalContent = 0; fNumberOfBatch = NumOfBatch; fKThreshold = CriticalityThreshold ; - SetBurnUpPrecision(0.005);//1 % of the targeted burnup + SetBurnUpPrecision(0.008);//0.8 % of the targeted burnup SetPCMPrecision(10); fInformationFile = InformationFile; @@ -160,7 +145,6 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(string TMVAWeightPath, int NumOfBatch, string Inform INFO << "\tThis model is based on the prediction of kinf" << endl; INFO << "\tThe TMVA (weight | information) files are :" << endl; INFO << "\t" << "( " << fTMVAWeightPath[0] << " | " << fInformationFile << " )" << endl; - INFO << "Maximal fissile content (molar proportion) : " << fMaximalContent << endl; INFO << "Maximal burnup (GWd/tHM) : " << fMaximalBU << endl; EquivalenceModel::PrintInfo(); @@ -174,13 +158,7 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(string TMVAWeightPath, int NumOfBatch, string Inform if(fStreamList[(* it_s_IV).first].GetIsotopicQuantity().empty()){EmptyList=true;} } - for( it_s_D = fFirstGuessContent .begin(); it_s_D != fFirstGuessContent.end(); it_s_D++) - { - if(fFirstGuessContent[(* it_s_D).first] ==0){FirstContentNULL=true;} - } - - - if(fMapOfTMVAVariableNames.empty() || EmptyList==true || FirstContentNULL==true || fMaximalBU == 0 || fMaximalContent == 0 ) + if(fMapOfTMVAVariableNames.empty() || EmptyList==true || fMaximalBU == 0 ) { ERROR<<"Missing information file in : "<<fInformationFile<<endl; exit(1); @@ -198,14 +176,14 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(CLASSLogger* log, string TMVAWeightPath, int NumOfBa InformationFile = StringLine::ReplaceAll(TMVAWeightPath,".xml",".nfo"); fMaximalBU = 0; - fMaximalContent = 0; fInformationFile = InformationFile; LoadKeyword(); ReadNFO();//Getting information from fMLPInformationFile fNumberOfBatch = NumOfBatch; fKThreshold = CriticalityThreshold ; - SetBurnUpPrecision(0.005);//1 % of the targeted burnup + + SetBurnUpPrecision(0.008);//0.8 % of the targeted burnup SetPCMPrecision(10); INFO << "__An equivalence model has been define__" << endl; @@ -227,13 +205,8 @@ EQM_MLP_Kinf::EQM_MLP_Kinf(CLASSLogger* log, string TMVAWeightPath, int NumOfBa if(fStreamList[(* it_s_IV).first].GetIsotopicQuantity().empty()){EmptyList=true;} } - for( it_s_D = fFirstGuessContent .begin(); it_s_D != fFirstGuessContent.end(); it_s_D++) - { - if(fFirstGuessContent[(* it_s_D).first] ==0){FirstContentNULL=true;} - } - - if(fMapOfTMVAVariableNames.empty() || EmptyList==true || FirstContentNULL==true || fMaximalBU == 0 || fMaximalContent == 0 ) + if(fMapOfTMVAVariableNames.empty() || EmptyList==true || fMaximalBU == 0 ) { ERROR<<"Missing information file in : "<<fInformationFile<<endl; exit(1); @@ -291,22 +264,6 @@ void EQM_MLP_Kinf::ReadMaxBurnUp(const string &line) DBGL } //________________________________________________________________________ -void EQM_MLP_Kinf::ReadMaxFisContent(const string &line) -{ - DBGL - int pos = 0; - string keyword = tlc(StringLine::NextWord(line, pos, ' ')); - if( keyword != "k_maxfiscontent" ) // Check the keyword - { - ERROR << " Bad keyword : \"k_maxfiscontent\" not found !" << endl; - exit(1); - } - - fMaximalContent = atof(StringLine::NextWord(line, pos, ' ').c_str()); - - DBGL -} -//________________________________________________________________________ void EQM_MLP_Kinf::ReadLine(string line) { DBGL diff --git a/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.cxx b/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.cxx index d4ee6b5661238b8c5802e125be40c05a989a2ddd..7a51fd6c4c466817dca16cd5b14bbb925422fb79 100644 --- a/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.cxx +++ b/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.cxx @@ -31,14 +31,12 @@ EQM_MLP_PWR_MOxEUS::EQM_MLP_PWR_MOxEUS(string TMVAWeightPath, int NumOfBatch, double CriticalityThreshold, double MaximalPuMassContent):EquivalenceModel(new CLASSLogger("EQM_MLP_PWR_MOxEUS.log")) { - fNumberOfBatch = NumOfBatch; + fNumberOfBatch = NumOfBatch; fKThreshold = CriticalityThreshold ; fMaximalPuMassContent = MaximalPuMassContent; fTMVAWeightPath.push_back(TMVAWeightPath); - SetMinimalU5Enrichment(0.003); - ZAI U8(92,238,0); ZAI U5(92,235,0); ZAI Pu8(94,238,0); @@ -51,16 +49,16 @@ EQM_MLP_PWR_MOxEUS::EQM_MLP_PWR_MOxEUS(string TMVAWeightPath, int NumOfBatch, do fStreamList["EnrichmentList"] = U5*1; fStreamList["FertileList"] = U8*1; - fFirstGuessContent["PuList"] = 0.02; - fFirstGuessContent["EnrichmentList"] = (1-fFirstGuessContent["PuList"])*GetMinimalU5Enrichment(); - fFirstGuessContent["FertileList"] = 1 - fFirstGuessContent["PuList"] - fFirstGuessContent["EnrichmentList"]; + fStreamListEqMMassFractionMin["PuList"] = 0.0; + fStreamListEqMMassFractionMin["EnrichmentList"] = 0.025; + + fStreamListEqMMassFractionMax["PuList"] = 0.16; + fStreamListEqMMassFractionMax["EnrichmentList"] = 0.05; fSpecificPower = 34.24; - fMaximalPuMolarContent = 0.50; - fMaximalBU = 100; + fMaximalBU = 75; SetRelativMassPrecision(0.0008); - SetMaximalU5Enrichment(0.10); SetBurnUpPrecision(0.008);//1 % of the targeted burnup SetPCMPrecision(10); @@ -80,8 +78,6 @@ EQM_MLP_PWR_MOxEUS::EQM_MLP_PWR_MOxEUS(CLASSLogger* log, string TMVAWeightPath, fTMVAWeightPath.push_back(TMVAWeightPath); - SetMinimalU5Enrichment(0.003); - ZAI U8(92,238,0); ZAI U5(92,235,0); ZAI Pu8(94,238,0); @@ -94,16 +90,16 @@ EQM_MLP_PWR_MOxEUS::EQM_MLP_PWR_MOxEUS(CLASSLogger* log, string TMVAWeightPath, fStreamList["EnrichmentList"] = U5*1; fStreamList["FertileList"] = U8*1; - fFirstGuessContent["PuList"] = 0.02; - fFirstGuessContent["EnrichmentList"] = (1-fFirstGuessContent["PuList"])*GetMinimalU5Enrichment(); - fFirstGuessContent["FertileList"] = 1 - fFirstGuessContent["PuList"] - fFirstGuessContent["EnrichmentList"]; + fStreamListEqMMassFractionMin["PuList"] = 0.0; + fStreamListEqMMassFractionMin["EnrichmentList"] = 0.025; + + fStreamListEqMMassFractionMax["PuList"] = 0.16; + fStreamListEqMMassFractionMax["EnrichmentList"] = 0.05; fSpecificPower = 34.24; - fMaximalPuMolarContent = 1.0; - fMaximalBU = 100; + fMaximalBU = 75; SetRelativMassPrecision(0.0008); - SetMaximalU5Enrichment(0.10); SetBurnUpPrecision(0.008);//1 % of the targeted burnup SetPCMPrecision(10); @@ -111,276 +107,6 @@ EQM_MLP_PWR_MOxEUS::EQM_MLP_PWR_MOxEUS(CLASSLogger* log, string TMVAWeightPath, INFO<<"\tThis model is based on a multi layer perceptron"<<endl; INFO<<"\t\tThe TMVA weight file is :"<<endl; INFO<<"\t\t\t"<<fTMVAWeightPath[0]<<endl; - -} - -//________________________________________________________________________ -map <string , vector<double> > EQM_MLP_PWR_MOxEUS::BuildFuel(double BurnUp, double HMMass, map < string , vector <IsotopicVector> > StreamArray) -{ - - map <string , vector<double> > lambda ; // map containing name of the list and associated vector of proportions taken from stocks - - //Iterators declaration - map < string , vector <IsotopicVector> >::iterator it_s_vIV; - map < string , vector <double> >::iterator it_s_vD; - map < string , IsotopicVector >::iterator it_s_IV; - map < string , double >::iterator it_s_D; - map < string , bool >::iterator it_s_B; - - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - { - for( size_t i=0; i<StreamArray[(*it_s_vIV).first].size(); i++) - { - lambda[(*it_s_vIV).first].push_back(0); - } - } - - /*** Test if there is stocks **/ - bool BreakReturnLambda = false; - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - { - if(StreamArray[(*it_s_vIV).first].size() == 0) - { - WARNING << " No stock available for stream : "<< (*it_s_vIV).first <<". Fuel not build." << endl; - SetLambdaToErrorCode(lambda[(*it_s_vIV).first]); - BreakReturnLambda = true; - } - } - if(BreakReturnLambda) { return lambda;} - HMMass *= 1e6; //Unit onversion : tons to gram - - /**** Some initializations **/ - - StocksTotalMassCalculation(StreamArray); - - fActualMassContentInFuel = GetBuildFuelFirstGuess(); - int loopCount = 0; - bool FuelBuiltCorrectly = false ; - - map <string , IsotopicVector > IVStream; - map <string , double > MaterialMolarContent; - map <string , double > MaterialMassContent; - map <string , double > MaterialMassNeeded; - map <string , double > LambdaNeeded; - map <string , double > DeltaMass; - - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - { - LambdaNeeded[(*it_s_vIV).first] = 0; - DeltaMass[(*it_s_vIV).first] = 0; - } - - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - { - MaterialMassNeeded[(*it_s_vIV).first] = HMMass*fActualMassContentInFuel[(*it_s_vIV).first]; - DeltaMass[(*it_s_vIV).first] = - MaterialMassNeeded[(*it_s_vIV).first]; - } - - do - { - map < string , double > MaterialMeanMolar; - map < string , double > MaterialMassAvailable; - map < string , bool > CheckOnMass; - - map < string , double > LambdaPreviousStep; - - double FuelMeanMolar = 0; - - bool NotEnoughPuInStock = false; - BreakReturnLambda = false; - - for( it_s_D = LambdaNeeded.begin(); it_s_D != LambdaNeeded.end(); it_s_D++) - LambdaPreviousStep[(*it_s_D ).first] =LambdaNeeded[(*it_s_D ).first]; - - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - LambdaNeeded[(*it_s_vIV).first] = LambdaCalculation((*it_s_vIV).first, LambdaPreviousStep[(*it_s_vIV).first], MaterialMassNeeded[(*it_s_vIV).first], DeltaMass[(*it_s_vIV).first], StreamArray[(*it_s_vIV).first]); - - - for( it_s_D = LambdaNeeded.begin(); it_s_D != LambdaNeeded.end(); it_s_D++) - { - if( LambdaNeeded[(*it_s_D).first] == -1 ) - { - SetLambdaToErrorCode(lambda[(*it_s_D).first]); - WARNING << "Not enough : "<< (*it_s_D).first <<" to build fuel." << endl; - BreakReturnLambda = true; - if((*it_s_D).first=="PuList") - NotEnoughPuInStock = true; - } - } - - if(BreakReturnLambda && !NotEnoughPuInStock){ return lambda;} - - if(BreakReturnLambda && NotEnoughPuInStock) - { - MaterialMolarContent["PuList"] = -1; - break; - } - - for( it_s_D = LambdaNeeded.begin(); it_s_D != LambdaNeeded.end(); it_s_D++) - SetLambda(lambda[(*it_s_D).first], LambdaNeeded[(*it_s_D).first]); - - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - IVStream[(*it_s_vIV).first].Clear(); - - for( it_s_vIV = StreamArray.begin(); it_s_vIV != StreamArray.end(); it_s_vIV++) - { - for(size_t i=0; i<StreamArray[(*it_s_vIV).first].size(); i++) - IVStream[(*it_s_vIV).first] += lambda[(*it_s_vIV).first][i] * StreamArray[(*it_s_vIV).first][i]; - } - - if (loopCount > fMaxInterration) - { - ERROR << "Too much iterration in BuildFuel Method !"; - ERROR << "Need improvement in fuel fabrication ! Ask for it or D.I.Y. !!" << endl; - exit(1); - } - - /* Calcul the quantity of this composition needed to reach the burnup */ - MaterialMolarContent = GetMolarFraction(IVStream, BurnUp); - //cout<<MaterialMolarContent["PuList"]<<endl; - for( it_s_D = MaterialMolarContent.begin(); it_s_D != MaterialMolarContent.end(); it_s_D++) - { - if( MaterialMolarContent[(*it_s_D).first] < 0 || MaterialMolarContent[(*it_s_D).first] > 1 ) - { - SetLambdaToErrorCode(lambda[(*it_s_D).first]); - WARNING << "GetMolarFraction return negative or greater than one value, at least for one material : "<<(*it_s_D).first; - return lambda; - } - } - for( it_s_IV = IVStream.begin(); it_s_IV != IVStream.end(); it_s_IV++) - MaterialMassAvailable[(*it_s_IV).first] = IVStream[(*it_s_IV).first].GetTotalMass()*1e06; - - for( it_s_D = MaterialMolarContent.begin(); it_s_D != MaterialMolarContent.end(); it_s_D++) - { - MaterialMeanMolar[(*it_s_D).first] = IVStream[(*it_s_D).first].GetMeanMolarMass(); - FuelMeanMolar += MaterialMolarContent[(*it_s_D).first] * MaterialMeanMolar[(*it_s_D).first]; - } - - for( it_s_D = MaterialMolarContent.begin(); it_s_D != MaterialMolarContent.end(); it_s_D++) - MaterialMassContent [(*it_s_D).first] = MaterialMolarContent[(*it_s_D).first] * MaterialMeanMolar[(*it_s_D).first] / FuelMeanMolar; - - fActualMolarContentInFuel = MaterialMolarContent; //fActualMolarContentInFuel is used in GetMolarFraction - - for( it_s_D = MaterialMassContent.begin(); it_s_D != MaterialMassContent.end(); it_s_D++) - { - MaterialMassNeeded[(*it_s_D).first] = HMMass*MaterialMassContent[(*it_s_D).first]; - DeltaMass[(*it_s_D).first] = MaterialMassAvailable[(*it_s_D).first] - MaterialMassNeeded[(*it_s_D).first]; - } - - for( it_s_D = MaterialMassNeeded.begin(); it_s_D != MaterialMassNeeded.end(); it_s_D++) - { - double DeltaM = fabs(MaterialMassNeeded[(*it_s_D).first] - MaterialMassAvailable[(*it_s_D).first]) / HMMass ; - if(DeltaM<fRelativMassPrecision) {CheckOnMass[(*it_s_D).first] = true;} - else{CheckOnMass[(*it_s_D).first] = false;} - } - - FuelBuiltCorrectly = true; - for( it_s_B = CheckOnMass.begin(); it_s_B != CheckOnMass.end(); it_s_B++) - { - if(!CheckOnMass[(*it_s_B).first]) {FuelBuiltCorrectly = false;} - } - - loopCount++; - - }while(!FuelBuiltCorrectly); - - if(MaterialMassContent["PuList"]>fMaximalPuMassContent || MaterialMolarContent["PuList"]==-1) - { - for( it_s_vD = lambda.begin(); it_s_vD != lambda.end(); it_s_vD++) - { - for(size_t i=0; i<lambda[(*it_s_vD).first].size(); i++) - { - lambda[(*it_s_vD).first][i] = 0; - } - } - - for( it_s_D = LambdaNeeded.begin(); it_s_D != LambdaNeeded.end(); it_s_D++) - LambdaNeeded[(*it_s_D).first] = 0; - - for( it_s_IV = IVStream.begin(); it_s_IV != IVStream.end(); it_s_IV++) - IVStream[(*it_s_IV).first] = IsotopicVector(); - - if(MaterialMassContent["PuList"]>fMaximalPuMassContent) - { - MaterialMolarContent["PuList"] = 0; - MaterialMassNeeded["PuList"] = fMaximalPuMassContent*HMMass; - LambdaNeeded["PuList"] = LambdaCalculation("PuList", 0., MaterialMassNeeded["PuList"], 1., StreamArray["PuList"]); - } - - if(MaterialMolarContent["PuList"]==-1) - { - LambdaNeeded["PuList"] = LambdaCalculation("PuList", 0., fTotalMassInStocks["PuList"], 1., StreamArray["PuList"]); - } - - SetLambda(lambda["PuList"], LambdaNeeded["PuList"]); - - for( int i = 0 ; i < (int)lambda["PuList"].size() ; i++ ) - IVStream["PuList"] +=lambda["PuList"][i] * StreamArray["PuList"][i]; - - fActualMassContentInFuel["PuList"] = (IVStream["PuList"].GetTotalMass()*1e6)/HMMass; - - LambdaNeeded["EnrichmentList"] = LambdaCalculation("EnrichmentList", 0., HMMass, 1., StreamArray["EnrichmentList"]); - LambdaNeeded["FertileList"] = LambdaCalculation("FertileList", 0., HMMass, 1., StreamArray["FertileList"]); - - SetLambda(lambda["EnrichmentList"], LambdaNeeded["EnrichmentList"]); - SetLambda(lambda["FertileList"], LambdaNeeded["FertileList"]); - - for( int i = 0 ; i < (int)lambda["EnrichmentList"].size() ; i++ ) - IVStream["EnrichmentList"] +=lambda["EnrichmentList"][i] * StreamArray["EnrichmentList"][i]; - - for( int i = 0 ; i < (int)lambda["FertileList"].size() ; i++ ) - IVStream["FertileList"] +=lambda["FertileList"][i] * StreamArray["FertileList"][i]; - - double U5Enrichment = GetU5Enrichment(IVStream, BurnUp); - double MeanMolarUdepleted = U5Enrichment*IVStream["EnrichmentList"].GetMeanMolarMass() + (1 - U5Enrichment)*IVStream["FertileList"].GetMeanMolarMass(); - - double UdepletedMassNeeded = (HMMass - IVStream["PuList"].GetTotalMass() * 1e6); - double EnrichedSupportMassNeeded = U5Enrichment*UdepletedMassNeeded*IVStream["EnrichmentList"].GetMeanMolarMass()/IVStream["FertileList"].GetMeanMolarMass(); - double FertilMassNeeded = UdepletedMassNeeded - EnrichedSupportMassNeeded; - - LambdaNeeded["EnrichmentList"] = LambdaCalculation("EnrichmentList", 0., EnrichedSupportMassNeeded, 1., StreamArray["EnrichmentList"]); - LambdaNeeded["FertileList"] = LambdaCalculation("FertileList", 0., FertilMassNeeded, 1., StreamArray["FertileList"]); - - BreakReturnLambda = false; - - for( it_s_D = LambdaNeeded.begin(); it_s_D != LambdaNeeded.end(); it_s_D++) - { - if( LambdaNeeded[(*it_s_D).first] == -1 ) - { - SetLambdaToErrorCode(lambda[(*it_s_D).first]); - WARNING << "Not enough : "<< (*it_s_D).first <<" to build fuel." << endl; - BreakReturnLambda = true; - } - } - if(BreakReturnLambda) { return lambda;} - - SetLambda(lambda["EnrichmentList"], LambdaNeeded["EnrichmentList"]); - SetLambda(lambda["FertileList"], LambdaNeeded["FertileList"]); - - IVStream["EnrichmentList"] = IsotopicVector(); - IVStream["FertileList"] = IsotopicVector(); - - for( int i = 0 ; i < (int)lambda["EnrichmentList"].size() ; i++ ) - IVStream["EnrichmentList"] +=lambda["EnrichmentList"][i] * StreamArray["EnrichmentList"][i]; - - for( int i = 0 ; i < (int)lambda["FertileList"].size() ; i++ ) - IVStream["FertileList"] +=lambda["FertileList"][i] * StreamArray["FertileList"][i]; - } - - for( it_s_D = MaterialMassNeeded.begin(); it_s_D != MaterialMassNeeded.end(); it_s_D++) - DBGV( "Weight percent : "<<(*it_s_D).first<<" "<< MaterialMassNeeded[(*it_s_D).first]/HMMass); - - for( it_s_vD = lambda.begin(); it_s_vD != lambda.end(); it_s_vD++) - { - DBGV( "Lambda vector : "<<(*it_s_vD).first ); - - for(size_t i=0; i<lambda[(*it_s_vD).first].size(); i++) - { - DBGV(lambda[(*it_s_vD).first][i]); - } - } - - return lambda; } //________________________________________________________________________ @@ -556,124 +282,3 @@ double EQM_MLP_PWR_MOxEUS::GetMaximumBurnUp(IsotopicVector TheFuel, double Targe return SecondToBurnup(TheFinalTime); } -//________________________________________________________________________ -map < string , double> EQM_MLP_PWR_MOxEUS::GetMolarFraction(map <string , IsotopicVector > IVStream, double TargetBU) -{ - - IsotopicVector Pu = IVStream["PuList"]; - IsotopicVector EnrichedSupport = IVStream["EnrichmentList"]; - IsotopicVector Fertil = IVStream["FertileList"]; - - map < string , double> MolarFraction ; - - double OldPuContentMinus = 0; - double OldPuContentPlus = fMaximalPuMolarContent; - - double PredictedBU = 0 ; - double PuContent = fActualMolarContentInFuel["PuList"]; - double U5Enrichment = fMinimalU5Enrichment; - - IsotopicVector FreshFuel = PuContent*(Pu/Pu.GetSumOfAll()) + U5Enrichment*(1-PuContent) *(EnrichedSupport/EnrichedSupport.GetSumOfAll()) + (1 - U5Enrichment)*(1-PuContent)*(Fertil/Fertil.GetSumOfAll()); - double OldPredictedBU = GetMaximumBurnUp(FreshFuel,TargetBU); - - double Precision = GetBurnUpPrecision()*TargetBU; //1 % of the targeted burnup - int count = 0; - int MaximumLoopCount = 500; - do - { - if(count > MaximumLoopCount ) - { - ERROR<<"CRITICAL ! Can't manage to predict fissile content\nHint : Try to decrease the precision on burnup using :\nYourEQM_MLP_Kinf->SetBurnUpPrecision(prop); with prop the precision (default 0.5percent : 0.005) INCREASE IT\nIf this message still appear mail to leniau@subatech.in2p3.fr\nor nicolas.thiolliere@subatech.in2p3.fr "<<endl; - exit(1); - } - if( (OldPredictedBU - TargetBU) < 0 ) //The Content can be increased - { - OldPuContentMinus = PuContent; - PuContent = PuContent + fabs(OldPuContentPlus-PuContent)/2.; - } - else if( (OldPredictedBU - TargetBU) > 0) //The Content is too high - { - OldPuContentPlus = PuContent; - PuContent = PuContent - fabs(OldPuContentMinus-PuContent)/2.; - } - - IsotopicVector FreshFuel= PuContent*(Pu/Pu.GetSumOfAll()) + U5Enrichment*(1-PuContent) *(EnrichedSupport/EnrichedSupport.GetSumOfAll()) + (1 - U5Enrichment)*(1-PuContent)*(Fertil/Fertil.GetSumOfAll()); - PredictedBU = GetMaximumBurnUp(FreshFuel,TargetBU); - - OldPredictedBU = PredictedBU; - count ++; -INFO<<"count "<<count<<" => Target BU "<<TargetBU<<" - Predicted BU "<<PredictedBU<<endl; - - }while(fabs(TargetBU-PredictedBU)>Precision); - - DBGV("Predicted BU "<<PredictedBU<<" PuContent "<<PuContent); - - MolarFraction["PuList"] = PuContent; - MolarFraction["EnrichmentList"] = (1 - MolarFraction["PuList"])*fMinimalU5Enrichment; - MolarFraction["FertileList"] = 1 - MolarFraction["PuList"] - MolarFraction["EnrichmentList"]; - -return MolarFraction; -} -//________________________________________________________________________ - -//________________________________________________________________________ -double EQM_MLP_PWR_MOxEUS::GetU5Enrichment(map <string , IsotopicVector > IVStream, double TargetBU) -{ - IsotopicVector Pu = IVStream["PuList"]; - IsotopicVector EnrichedSupport = IVStream["EnrichmentList"]; - IsotopicVector Fertil = IVStream["FertileList"]; - - double OldU5EnrichmentMinus = fMinimalU5Enrichment; - double OldU5EnrichmentPlus = fMaximalU5Enrichment; - double U5Enrichment = fMinimalU5Enrichment; - - double PredictedBU = 0 ; - double USupportMassNeeded = (Pu.GetTotalMass()*1e6)/fActualMassContentInFuel["PuList"] - (Pu.GetTotalMass()*1e6); - - IsotopicVector USupport = U5Enrichment*EnrichedSupport + (1 - U5Enrichment)*Fertil; - double USupportMassContent = USupportMassNeeded/(USupport.GetTotalMass()*1e6); - - IsotopicVector FreshFuel = Pu + USupportMassContent*USupport; - FreshFuel = FreshFuel/FreshFuel.GetSumOfAll(); - - double OldPredictedBU = GetMaximumBurnUp(FreshFuel,TargetBU); - - double Precision = GetBurnUpPrecision()*TargetBU; //1 % of the targeted burnup - int count = 0; - int MaximumLoopCount = 500; - - do - { - - if(count > MaximumLoopCount ) - { - ERROR<<"CRITICAL ! Can't manage to predict fissile content\nHint : Try to decrease the precision on burnup using :\nYourEQM_MLP_Kinf->SetBurnUpPrecision(prop); with prop the precision (default 0.5percent : 0.005) INCREASE IT\nIf this message still appear mail to leniau@subatech.in2p3.fr\nor nicolas.thiolliere@subatech.in2p3.fr "<<endl; - exit(1); - } - if( (OldPredictedBU - TargetBU) < 0 ) //The Content can be increased - { - OldU5EnrichmentMinus = U5Enrichment; - U5Enrichment = U5Enrichment + fabs(OldU5EnrichmentPlus-U5Enrichment)/2.; - } - else if( (OldPredictedBU - TargetBU) > 0) //The Content is too high - { - OldU5EnrichmentPlus = U5Enrichment; - U5Enrichment = U5Enrichment - fabs(OldU5EnrichmentMinus-U5Enrichment)/2.; - } - - IsotopicVector USupport = U5Enrichment*EnrichedSupport + (1 - U5Enrichment)*Fertil; - - double USupportMassContent = USupportMassNeeded/(USupport.GetTotalMass()*1e6); - IsotopicVector FreshFuel = Pu + USupportMassContent*USupport; - FreshFuel = FreshFuel/FreshFuel.GetSumOfAll(); - - PredictedBU = GetMaximumBurnUp(FreshFuel,TargetBU); - - OldPredictedBU = PredictedBU; - count ++; - }while(fabs(TargetBU-PredictedBU)>Precision); - -DBGV("Predicted BU "<<PredictedBU<<" U5Enrichment "<<U5Enrichment); -return U5Enrichment; - -} diff --git a/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.hxx b/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.hxx index e4a138cf3bdc28efd671e62954162703a36a8ee6..96aaa0d1fe1fcac3fdc9685c9a89ffa7bd28c5c5 100644 --- a/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.hxx +++ b/source/Model/Equivalence/EQM_MLP_PWR_MOxEUS.hxx @@ -22,21 +22,14 @@ class EQM_MLP_PWR_MOxEUS : public EquivalenceModel virtual map <string , vector<double> > BuildFuel(double BurnUp, double HMMass, map < string , vector <IsotopicVector> > StreamArray); void SetSpecificPower(double SpecificPower) {fSpecificPower = SpecificPower;} - void SetMinimalU5Enrichment(double MinEU5) {fMinimalU5Enrichment= MinEU5;} - void SetMaximalU5Enrichment(double MaxEU5){fMaximalU5Enrichment= MaxEU5;} - void SetBurnUpPrecision(double precision) {fBurnUpPrecision= precision;} void SetPCMPrecision(double prop) {fPCMPrecision= prop;} - double GetBurnUpPrecision(){return fBurnUpPrecision;} double GetPCMPrecision(){return fPCMPrecision/1e5;} double GetMinimalU5Enrichment(){return fMinimalU5Enrichment;}; TTree* CreateTMVAInputTree(IsotopicVector TheFuel, double ThisTime); double ExecuteTMVA(TTree* theTree, string WeightPath); double GetMaximumBurnUp(IsotopicVector TheFuel, double TargetBU); - double GetU5Enrichment(map <string , IsotopicVector > IVStream, double TargetBU); - - map < string , double> GetMolarFraction(map < string , IsotopicVector> IVStream, double BurnUp); private: @@ -47,15 +40,7 @@ class EQM_MLP_PWR_MOxEUS : public EquivalenceModel double fKThreshold ; double fSpecificPower ; double fMaximalBU; - double fBurnUpPrecision; double fPCMPrecision; - - double fActualPuMolarContent; - double fActualPuMassContent; - double fMaximalPuMassContent ; - double fMaximalPuMolarContent; - double fMinimalU5Enrichment; - double fMaximalU5Enrichment; vector <string> fTMVAWeightPath; //!<The weight needed by TMVA to construct and execute the multilayer perceptron diff --git a/source/include/EquivalenceModel.hxx b/source/include/EquivalenceModel.hxx index a5cc8eb7b19d15d50e5088fd0631a18e7afe5cef..82d560fb2f080d4efb87746d0fe868fcc549aaa3 100644 --- a/source/include/EquivalenceModel.hxx +++ b/source/include/EquivalenceModel.hxx @@ -94,11 +94,11 @@ class EquivalenceModel : public CLASSObject map < string, IsotopicVector> GetAllStreamList() {return fStreamList;} //!<return all the lists int GetStreamListNumber(){return fStreamList.size();}; - int GetMaxInterration() const { return fMaxInterration; } //!< Max iterration in build fueld algorythm + int GetMaxIterration() const { return fMaxInterration; } //!< Max iterration in build fueld algorythm double GetBurnUpPrecision(){return fBurnUpPrecision;}//!< Get the precision on Burnup : proportion of the targeted burnup - void SetRelativMassPrecision( double val) { fRelativMassPrecision = val; } //!< Mass precision - void SetMaxInterration(int val) { fMaxInterration = val; } //!< Max iteration in build fuel algorithm + void SetMaxIterration(int val) { fMaxInterration = val; } //!< Max iteration in build fuel algorithm + void SetBurnUpPrecision(double prop){fBurnUpPrecision = prop;} //!< Set the precision on Burnup : proportion of the targeted burnup //@} @@ -141,7 +141,8 @@ class EquivalenceModel : public CLASSObject */ void ReadList(const string &line); - void ReadFirstGuessContent(const string &line); + void ReadEqMaxFraction(const string &line); + void ReadEqMinFraction(const string &line); bool isIVInDomain(IsotopicVector IV); @@ -154,20 +155,14 @@ class EquivalenceModel : public CLASSObject //!< each list is identified by a keyword (example : -> "Fissile" & "Fertile") map < string , double> fStreamListEqMMassFractionMax; //!< Map that contains lists of stream according to the EqModel with mass maximum fraction map < string , double> fStreamListEqMMassFractionMin; //!< Map that contains lists of stream according to the EqModel with mass minimum fraction - - map < string, double> fFirstGuessContent; //!< fissile content for BuildFuel initialization (in weight proportion) - map < string, double> fActualMolarContentInFuel; //!< Molar Content in fuel of each list at this step of the calculation - map < string, double> fActualMassContentInFuel; //!< Mass Content in fuel of each list at this step of the calculation - double fSpecificPower; //!< The specific power in W/gHM (HM: heavy Metal) - double fMaximalContent; //!< The approx. maximum fissile content reachable by the model + double fSpecificPower; //!< The specific power in W/gHM (HM: heavy Metal) - void SetLambdaToErrorCode(vector<double>& lambda); - double fRelativMassPrecision; //!< Mass precision double fBurnUpPrecision; //!< precision on Burnup int fMaxIterration; //!< Max iterration in build fueld algorythm + void SetLambdaToErrorCode(vector<double>& lambda); //@} @@ -182,17 +177,6 @@ class EquivalenceModel : public CLASSObject string fDBFType; //!< Fuel Type (e.g MOX, UOX, ThU, ThPu ...) string fDBRType; //!< Reactor Type (e.g PWR, FBR-Na, ADS..) - /*! - \name Others - */ - //@{ - map <string , double > fTotalMassInStocks; //!< Total mass in each vector of stock - map <string , double > fLambdaMax; //!< Total lambda of available stocks - - - - - //@} }; #endif diff --git a/source/src/EquivalenceModel.cxx b/source/src/EquivalenceModel.cxx index 28310e69c7dd147e26c7f14a40bec78844eb7224..c44eb82b16449a45fc5f13e0bc212c909bcca1cc 100644 --- a/source/src/EquivalenceModel.cxx +++ b/source/src/EquivalenceModel.cxx @@ -5,7 +5,6 @@ //________________________________________________________________________ EquivalenceModel::EquivalenceModel():CLASSObject() { - fRelativMassPrecision = 5/10000.; //Mass precision fMaxIterration = 500; // Max iterration in build fueld algorythum freaded = false; @@ -14,7 +13,6 @@ EquivalenceModel::EquivalenceModel():CLASSObject() //________________________________________________________________________ EquivalenceModel::EquivalenceModel(CLASSLogger* log):CLASSObject(log) { - fRelativMassPrecision = 5/10000.; //Mass precision fMaxIterration = 500; // Max iterration in build fueld algorythm freaded = false; @@ -188,22 +186,41 @@ void EquivalenceModel::ReadList(const string &line) DBGL } //________________________________________________________________________ -void EquivalenceModel::ReadFirstGuessContent(const string &line) +void EquivalenceModel::ReadEqMinFraction(const string &line) { DBGL int pos = 0; string keyword = tlc(StringLine::NextWord(line, pos, ' ')); - if( keyword != "k_firstguesscontent" ) // Check the keyword + if( keyword != "k_MassFractionMin" ) // Check the keyword { - ERROR << " Bad keyword : \"k_firstguesscontent\" not found !" << endl; + ERROR << " Bad keyword : \"k_MassFractionMin\" not found !" << endl; exit(1); } - string ListName = StringLine::NextWord(line, pos, ' '); - double Q = atof(StringLine::NextWord(line, pos, ' ').c_str()); - fFirstGuessContent[ListName] = Q; - + string ListName= StringLine::NextWord(line, pos, ' '); + double Q = atof(StringLine::NextWord(line, pos, ' ').c_str()); + fStreamListEqMMassFractionMin[ListName] = Q; + DBGL } + +//________________________________________________________________________ +void EquivalenceModel::ReadEqMaxFraction(const string &line) +{ + DBGL + int pos = 0; + string keyword = tlc(StringLine::NextWord(line, pos, ' ')); + if( keyword != "k_MassFractionMax" ) // Check the keyword + { + ERROR << " Bad keyword : \"k_MassFractionMax\" not found !" << endl; + exit(1); + } + string ListName= StringLine::NextWord(line, pos, ' '); + double Q = atof(StringLine::NextWord(line, pos, ' ').c_str()); + fStreamListEqMMassFractionMax[ListName] = Q; + + DBGL +} + //________________________________________________________________________ void EquivalenceModel::ReadSpecificPower(const string &line) {