From 851036fc724520cb650fda544c8cb05a40d222ec Mon Sep 17 00:00:00 2001 From: Baptiste Mouginot <mouginot.baptiste@gmail.com> Date: Tue, 19 Mar 2013 10:17:18 +0000 Subject: [PATCH] CLASS Version 1.3 : CLASS, REACTOR, EVOLUTIONDATABASE, FABRICATIONPLANT : Add Warning. EVOLUTIONDATABASE : Debug evolution Modul (seams to be functional, need more extensive test) FABRICATIONPLANT : Add choice to Use fixed reference DB (for the new fuel Evolution Calculation) or the update at each evolution step (via the fUpdateReferenceDBatEachStep boolean ) REACTOR : Add GetPower method, and some comment (unit) EVOLUTIVEPRODUCT : ADD GetPower and GetHMMass method, and fCycleTime variable (and correct the reading of CycleTime) git-svn-id: svn+ssh://svn.in2p3.fr/class@65 0e7d625b-0364-4367-a6be-d5be4a48d228 --- source/trunk/include/CLASS.hxx | 46 ++-- source/trunk/include/EvolutiveProduct.hxx | 9 +- source/trunk/include/FabricationPlant.hxx | 3 +- source/trunk/include/Reactor.hxx | 5 +- source/trunk/src/CLASS.cxx | 47 +++- source/trunk/src/EvolutionDataBase.cxx | 316 +++++++++++++--------- source/trunk/src/EvolutiveProduct.cxx | 10 +- source/trunk/src/FabricationPlant.cxx | 154 +++++++---- source/trunk/src/Reactor.cxx | 86 ++++++ 9 files changed, 446 insertions(+), 230 deletions(-) diff --git a/source/trunk/include/CLASS.hxx b/source/trunk/include/CLASS.hxx index 483479d09..fc57e1080 100755 --- a/source/trunk/include/CLASS.hxx +++ b/source/trunk/include/CLASS.hxx @@ -112,7 +112,7 @@ public : protected : - LogFile* fLog; + LogFile* fLog; //!< Pointer to the Log bool fStockManagement; ///< True if real StockManagement false unstead @@ -128,28 +128,28 @@ protected : ///< 16 Fuel Fabrication - vector<Storage*> fStorage; ///< Vector of Storages - vector<TreatmentFactory*> fTreatmentFactory; ///< Vector of Treament Factory - vector<Reactor*> fReactor; ///< Vector of Reactor - vector<FabricationPlant*> fFabricationPlant; ///< Vector of FabricationPlant - EvolutionDataBase<ZAI>* fDecayDataBase; //!< Pointer to the Decay DataBase - - - TFile* fOutFile; ///< Pointer to the Root Output File - string fOutputFileName; //! Name of the Output File - TTree* fOutT; ///< Pointer to the Root Output TTree - string fOutputTreeName; //! Name of the Output TTree - - - IsotopicVector fWaste; ///< Waste IV - IsotopicVector fTotalStorage; ///< Sum of all IV in Storage IV - IsotopicVector fGod; ///< GodIncome IV - IsotopicVector fTotalCooling; ///< Sum of all IV in Cooling IV - IsotopicVector fFuelFabrication; ///< Sum of all IV in Fabrication IV - IsotopicVector fTotalInReactor; ///< Sum of all IV in Reactor IV - IsotopicVector fIVInCycleTotal; ///< Summ of all IV in the cycle (without Waste) IV - IsotopicVector fIVTotal; ///< Summ of all IV in the parc (including Waste) IV - double fParcPower; ///< Summ of the Power of all reactor in the parc + vector<Storage*> fStorage; ///< Vector of Storages + vector<TreatmentFactory*> fTreatmentFactory; ///< Vector of Treament Factory + vector<Reactor*> fReactor; ///< Vector of Reactor + vector<FabricationPlant*> fFabricationPlant; ///< Vector of FabricationPlant + EvolutionDataBase<ZAI>* fDecayDataBase; //!< Pointer to the Decay DataBase + + + TFile* fOutFile; ///< Pointer to the Root Output File + string fOutputFileName; //! Name of the Output File + TTree* fOutT; ///< Pointer to the Root Output TTree + string fOutputTreeName; //! Name of the Output TTree + string fOutputLogName; ///< Name of the Ouput log File + + IsotopicVector fWaste; ///< Waste IV + IsotopicVector fTotalStorage; ///< Sum of all IV in Storage IV + IsotopicVector fGod; ///< GodIncome IV + IsotopicVector fTotalCooling; ///< Sum of all IV in Cooling IV + IsotopicVector fFuelFabrication; ///< Sum of all IV in Fabrication IV + IsotopicVector fTotalInReactor; ///< Sum of all IV in Reactor IV + IsotopicVector fIVInCycleTotal; ///< Summ of all IV in the cycle (without Waste) IV + IsotopicVector fIVTotal; ///< Summ of all IV in the parc (including Waste) IV + double fParcPower; ///< Summ of the Power of all reactor in the parc }; diff --git a/source/trunk/include/EvolutiveProduct.hxx b/source/trunk/include/EvolutiveProduct.hxx index fa53dcd06..cbda6ffd4 100755 --- a/source/trunk/include/EvolutiveProduct.hxx +++ b/source/trunk/include/EvolutiveProduct.hxx @@ -69,7 +69,8 @@ public : TGraph* GetKeff() const { return fKeff; } TGraph* GetFlux() const { return fFlux; } - double GetPower() const { return fPower; } //!< + double GetPower() const { return fPower; } //!< + double GetHMMass() const { return fHMMass; } string GetDB_file() const { return fDB_file; } TGraph* GetEvolutionTGraph(const ZAI& zai); @@ -100,11 +101,13 @@ protected : - string fReactorType; + string fReactorType; string fFuelType; - double fPower; + double fPower; + double fCycleTime; double fHMMass; + void ReadDB(string DBfile); double Interpolate(double t, TGraph& EvolutionGraph); ///< Interpolating the value of EvolutionGraph at the t time diff --git a/source/trunk/include/FabricationPlant.hxx b/source/trunk/include/FabricationPlant.hxx index 3792c2e97..1155b853a 100644 --- a/source/trunk/include/FabricationPlant.hxx +++ b/source/trunk/include/FabricationPlant.hxx @@ -51,6 +51,7 @@ public : //********* Set Method *********// + void SetUpdateReferenceDBatEachStep(bool val){ fUpdateReferenceDBatEachStep = val;} void SetId(int id) { fId = id; } //!< Set The FB Parc'Id void SetParc(CLASS* parc) { fParc = parc; } //!< Set the Pointer to the Parc void SetLog(LogFile* Log) { fLog = Log; } //!< Set the Pointer to the Log @@ -105,7 +106,7 @@ public : protected : int fId; //!< Identity of the FabricationPlant inside the Parc cSecond fInternalTime; ///< Internal Clock - + bool fUpdateReferenceDBatEachStep; //********* Internal Parameter *********// CLASS* fParc; //!< Pointer to the Parc diff --git a/source/trunk/include/Reactor.hxx b/source/trunk/include/Reactor.hxx index 0d198d2f8..d055ccaf7 100755 --- a/source/trunk/include/Reactor.hxx +++ b/source/trunk/include/Reactor.hxx @@ -64,6 +64,7 @@ public : IsotopicVector GetIVBeginCycle() const { return fIVBeginCycle; } //!< Return the Starting Cycle IV (Note : IVBegin != IVIn, only if using charging plan) IsotopicVector GetIVOutCycle() const { return fIVOutCycle; } //!< Return the Out Cycle IV IsotopicVector GetIVInCycle() const { return fIVInCycle; } //!< Return the In Cycle IV (Note : IVIn != IVBegin, only if using charging plan) + double GetPower() const { return fPower; } //!< Return the cycle time of the Reactor cSecond GetCycleTime() const { return fCycleTime; } //!< Return the cycle time of the Reactor cSecond GetCreationTime() const { return fCreationTime; } //!< Return the creation time of the Reactor cSecond GetLifeTime() const { return fLifeTime; } //!< Return the creation time of the Reactor @@ -123,9 +124,9 @@ protected : EvolutionDataBase<IsotopicVector>* fFuelTypeDB; //! Pointer to a Fuel Type Database cSecond fCreationTime; ///< CLASS Universal Time of Creation - cSecond fLifeTime; ///< LifeTime Of the Reactor + cSecond fLifeTime; ///< LifeTime Of the Reactor (Operating's Duration) cSecond fCycleTime; ///< Cycle Time - double fPower; ///< Power + double fPower; ///< Power (in Watt) IsotopicVector fIVReactor; ///< Fuel evoluated IV in the reactor IsotopicVector fIVBeginCycle; ///< Fuel IV at the Beginning of a Cycle diff --git a/source/trunk/src/CLASS.cxx b/source/trunk/src/CLASS.cxx index 9c3873bfb..0092913af 100755 --- a/source/trunk/src/CLASS.cxx +++ b/source/trunk/src/CLASS.cxx @@ -50,14 +50,33 @@ CLASS::CLASS() DBGL; fPrintStep = (cSecond)(3600*24*365.25); // One Step per Year fAbsoluteTime = 0; - fStockManagement = true; fStartingTime = 0; + + fStockManagement = true; + fOutputFileName = "CLASS_Default.root"; fOutputTreeName = "Data"; - fLog = new LogFile("CLASS.log"); + + fOutputLogName = "CLASS.log"; + fLog = new LogFile(fOutputLogName.c_str()); fParcPower = 0; + // Warning + + cout << "!!Info!! !!!CLASS!!! A Parc has been define :" << endl; + cout << "\t Print set at : " << (double)(fPrintStep/3600/24/365.25) << " year" << endl; + cout << "\t Absolute Time set at " << (double)(fAbsoluteTime/3600/24/365.25) << " year" << endl; + cout << "\t StockManagement set at : true" << endl; + cout << "\t OutPut will be in \"" << fOutputFileName << "\" File and \"" << fOutputTreeName << "\" TTree" << endl; + cout << "\t Log will be in " << fOutputLogName << endl << endl; + + fLog->fLog << "!!Info!! !!!CLASS!!! Parc has been define :" << endl; + fLog->fLog << "\t Print set at : " << (double)(fPrintStep/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t StockManagement set at : true" << endl; + fLog->fLog << "\t OutPut will be in \"" << fOutputFileName << "\" File and \"" << fOutputTreeName << "\" TTree" << endl; + fLog->fLog << "\t Log will be in " << fOutputLogName << endl << endl; + DBGL; } @@ -68,13 +87,33 @@ CLASS::CLASS(double abstime) DBGL; fPrintStep = (cSecond)(3600*24*365.25); // One Step per Year fAbsoluteTime = (cSecond)abstime; - fStockManagement = true; fStartingTime = fAbsoluteTime; + + fStockManagement = true; + fOutputFileName = "CLASS_Default.root"; fOutputTreeName = "Data"; - fLog = new LogFile("CLASS.log"); + + fOutputLogName = "CLASS.log" ; + fLog = new LogFile(fOutputLogName.c_str()); fParcPower = 0; + + + + // Warning + + cout << "!!Info!! !!!CLASS!!! A Parc has been define :" << endl; + cout << "\t Print set at : " << (double)(fPrintStep/3600/24/365.25) << " year" << endl; + cout << "\t StockManagement set at : true" << endl; + cout << "\t OutPut will be in \"" << fOutputFileName << "\" File and \"" << fOutputTreeName << "\" TTree" << endl; + cout << "\t Log will be in " << fOutputLogName << endl; + + fLog->fLog << "!!Info!! !!!CLASS!!! Parc has been define :" << endl; + fLog->fLog << "\t Print set at : " << (double)(fPrintStep/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t StockManagement set at : true" << endl; + fLog->fLog << "\t OutPut will be in \"" << fOutputFileName << "\" File and \"" << fOutputTreeName << "\" TTree" << endl; + fLog->fLog << "\t Log will be in " << fOutputLogName << endl << endl; DBGL; diff --git a/source/trunk/src/EvolutionDataBase.cxx b/source/trunk/src/EvolutionDataBase.cxx index b2592acf7..db020622a 100755 --- a/source/trunk/src/EvolutionDataBase.cxx +++ b/source/trunk/src/EvolutionDataBase.cxx @@ -38,6 +38,16 @@ EvolutionDataBase<ZAI>::EvolutionDataBase(LogFile* Log, string DB_index_file) DBGL; fLog = Log; fDataBaseIndex = DB_index_file; + + + // Warning + + cout << "!!Info!! !!!EvolutionDataBase<ZAI>!!! A EvolutionData<ZAI> has been define :" << endl; + cout << "\t His index is : \"" << DB_index_file << "\"" << endl << endl; + + fLog->fLog << "!!Info!! !!!EvolutionDataBase<ZAI>!!! A EvolutionData<ZAI> has been define :" << endl; + fLog->fLog << "\t His index is : \"" << DB_index_file << "\"" << endl << endl; + DBGL; } @@ -54,9 +64,9 @@ IsotopicVector EvolutionDataBase<ZAI>::Evolution(const ZAI& zai, double dt) { DBGL; IsotopicVector returnIV; - + map<ZAI ,EvolutiveProduct >::iterator it = fEvolutionDataBase.find(zai); - + if (it == fEvolutionDataBase.end() ) { ifstream DB_index(fDataBaseIndex.c_str()); @@ -69,43 +79,43 @@ IsotopicVector EvolutionDataBase<ZAI>::Evolution(const ZAI& zai, double dt) bool zaifind = false; string tmp; getline(DB_index,tmp); // Read first line - while (!DB_index.eof()) + while (!DB_index.eof()) { string line; int start=0; getline(DB_index,line); // Read first line string first=StringLine::NextWord(line,start); // Read first word - + if(first.size()==0) break; // If First word is null.... quit - + int rZ=atoi(first.c_str()); // Get Z int rA=atoi(StringLine::NextWord(line,start).c_str()); // Get A int rI=atoi(StringLine::NextWord(line,start).c_str()); // Get Isomeric State - + if(rZ == zai.Z() && rA == zai.A() && rI == zai.I() ) { string file_name = StringLine::NextWord(line,start); EvolutiveProduct evolutionproduct = EvolutiveProduct(fLog, file_name); - #pragma omp critical(DBupdate) +#pragma omp critical(DBupdate) {fEvolutionDataBase.insert( pair<ZAI ,EvolutiveProduct >(zai, evolutionproduct) );} returnIV = evolutionproduct.GetIsotopicVectorAt(dt); zaifind = true; } } - - if(zaifind == false) + + if(zaifind == false) { { - fLog->fLog << "!!Warning!! !!!EVOLUTIVE DB!!! Oups... Can't Find the ZAI : " - << zai.Z() << " " << zai.A() << " " << zai.I() << "!!! It will be considered as stable !!" << endl; + fLog->fLog << "!!Warning!! !!!EVOLUTIVE DB!!! Oups... Can't Find the ZAI : " + << zai.Z() << " " << zai.A() << " " << zai.I() << "!!! It will be considered as stable !!" << endl; EvolutiveProduct evolutionproduct = EvolutiveProduct(fLog," " , true, zai); {fEvolutionDataBase.insert( pair<ZAI, EvolutiveProduct >(zai, evolutionproduct) );} returnIV = evolutionproduct.GetIsotopicVectorAt(dt); - + } } - - + + } else returnIV = (*it).second.GetIsotopicVectorAt(dt); @@ -117,11 +127,11 @@ bool EvolutionDataBase<ZAI>::IsDefine(const ZAI& zai) const { DBGL; map<ZAI ,EvolutiveProduct > evolutiondb = (*this).GetEvolutionDataBase(); - if (evolutiondb.find(zai) != evolutiondb.end()) + if (evolutiondb.find(zai) != evolutiondb.end()) return true; - else + else return false; - + } //________________________________________________________________________ //________________________________________________________________________ @@ -140,7 +150,19 @@ EvolutionDataBase<IsotopicVector>::EvolutionDataBase(LogFile* Log, string DB_ind DBGL; fLog = Log; fDataBaseIndex = DB_index_file; + ReadDataBase(); + + + // Warning + cout << "!!Info!! !!!EvolutionDataBase<IsotopicVector>!!! A EvolutionData<ZAI> has been define :" << endl; + cout << "\t His index is : \"" << DB_index_file << "\"" << endl; + cout << "\t " << fEvolutionDataBase.size() << " EvolutiveProduct have been read."<< endl << endl; + + fLog->fLog << "!!Info!! !!!EvolutionDataBase<IsotopicVector>!!! A EvolutionData<ZAI> has been define :" << endl; + fLog->fLog << "\t His index is : \"" << DB_index_file << "\"" << endl; + fLog->fLog << "\t " << fEvolutionDataBase.size() << " EvolutiveProduct have been read."<< endl << endl; + DBGL; } @@ -156,15 +178,15 @@ void EvolutionDataBase<IsotopicVector>::ReadDataBase() } vector<double> vTime; vector<double> vTimeErr; - + string line; int start = 0; // First Get Fuel Type - getline(DataDB, line); - if( StringLine::NextWord(line, start, ' ') != "TYPE") + getline(DataDB, line); + if( StringLine::NextWord(line, start, ' ') != "TYPE") { cout << "!!Bad Trouble!! !!!EvolutionDataBase!!! Bad Database file : " << fDataBaseIndex << " Can't find the type of the DataBase"<< endl; fLog->fLog << "!!Bad Trouble!! !!!EvolutionDataBase!!! Bad Database file : " << fDataBaseIndex << " Can't find the type of the DataBase"<< endl; @@ -172,9 +194,9 @@ void EvolutionDataBase<IsotopicVector>::ReadDataBase() } fFuelType = StringLine::NextWord(line, start, ' '); // First Get Fuel Parameter - getline(DataDB, line); + getline(DataDB, line); start = 0; - if( StringLine::NextWord(line, start, ' ') != "PARAM") + if( StringLine::NextWord(line, start, ' ') != "PARAM") { cout << "!!Bad Trouble!! !!!EvolutionDataBase!!! Bad Database file : " << fDataBaseIndex << " Can't find the Parameter of the DataBase"<< endl; fLog->fLog << "!!Bad Trouble!! !!!EvolutionDataBase!!! Bad Database file : " << fDataBaseIndex << " Can't find the Parameter of the DataBase"<< endl; @@ -182,13 +204,13 @@ void EvolutionDataBase<IsotopicVector>::ReadDataBase() } while(start < (int)line.size()) fFuelParameter.push_back(atof(StringLine::NextWord(line, start, ' ').c_str())); - + //Then Get All the Database - + while (!DataDB.eof()) { - getline(DataDB, line); + getline(DataDB, line); if(line != "") { EvolutiveProduct* evolutionproduct = new EvolutiveProduct(fLog, line); @@ -207,7 +229,7 @@ map<double, EvolutiveProduct> EvolutionDataBase<IsotopicVector>::GetDistancesTo( { DBGL; map<double, EvolutiveProduct> distances; - + map<IsotopicVector, EvolutiveProduct > evolutiondb = fEvolutionDataBase; map<IsotopicVector, EvolutiveProduct >::iterator it; @@ -215,10 +237,10 @@ map<double, EvolutiveProduct> EvolutionDataBase<IsotopicVector>::GetDistancesTo( { pair<map<double, EvolutiveProduct>::iterator, bool> IResult; double D = Distance(isotopicvector.GetActinidesComposition(), (*it).second.GetIsotopicVectorAt(t).GetActinidesComposition()/ Norme( (*it).second.GetIsotopicVectorAt(t).GetActinidesComposition() )*Norme(isotopicvector.GetActinidesComposition()) ); -// double D = RelativDistance(isotopicvector, (*it).second.GetIsotopicVectorAt(t)/ Norme( (*it).second.GetIsotopicVectorAt(t) )*Norme(isotopicvector) ); + // double D = RelativDistance(isotopicvector, (*it).second.GetIsotopicVectorAt(t)/ Norme( (*it).second.GetIsotopicVectorAt(t) )*Norme(isotopicvector) ); IResult = distances.insert( pair<double, EvolutiveProduct>( D , (*it).second ) ); } - + return distances; DBGL; } @@ -346,7 +368,7 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is map< ZAI, double > toAdd; toAdd.insert(pair<ZAI, double> ( ZAI(94,243,0) , 1) ); ZAIDecay.insert( pair< ZAI, pair<double, map< ZAI, double > > >( ZAI(96,247,0), pair<double, map< ZAI, double > > ( 4.92298560000000000e+14, toAdd) ) ); - } + } { // 248Cm map< ZAI, double > toAdd; toAdd.insert(pair<ZAI, double> ( ZAI(-3,-3,-3) , 1) ); @@ -358,23 +380,23 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(-3,-3,-3) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(90,231,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(90,231,0), toAdd ) ); } { // 233Th map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(92,233,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(90,233,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(90,233,0), toAdd ) ); } { // 233Pa map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(92,233,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(91,233,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(91,233,0), toAdd ) ); } { // 237U map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(93,237,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(92,237,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(92,237,0), toAdd ) ); } { // 239U map<ZAI, double> toAdd ; @@ -384,43 +406,43 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is { // 238Np map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(94,238,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,238,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,238,0), toAdd ) ); } { // 239Np map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(94,239,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,239,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,239,0), toAdd ) ); } { // 240Np map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(94,240,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,240,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,240,0), toAdd ) ); } { // 241Np map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(94,241,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,241,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(93,241,0), toAdd ) ); } { // 237Pu map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(93,237,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(94,237,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(94,237,0), toAdd ) ); } { // 243Pu map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(95,243,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(94,243,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(94,243,0), toAdd ) ); } { // 240Am map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(94,240,0) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,240,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,240,0), toAdd ) ); } { // 242Am map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(96,242,0) , 0.827) ); toAdd.insert(pair<ZAI, double> ( ZAI(94,242,0) , 0.173) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,242,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,242,0), toAdd ) ); } { // 244Am map<ZAI, double> toAdd ; @@ -435,21 +457,21 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is { // 249Cm map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(-3,-3,-3) , 1) ); - FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(96,249,0), toAdd ) ); + FastDecay.insert( pair< ZAI, map<ZAI, double> > ( ZAI(96,249,0), toAdd ) ); } - - + + map<ZAI, map<ZAI, double> > Capture; { // 241Am map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(95,242,0) , 0.086) ); toAdd.insert(pair<ZAI, double> ( ZAI(95,242,1) , 0.914) ); - Capture.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,241,0), toAdd ) ); + Capture.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,241,0), toAdd ) ); } { // 242Am* map<ZAI, double> toAdd ; toAdd.insert(pair<ZAI, double> ( ZAI(95,243,0) , 1) ); - Capture.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,242,1), toAdd ) ); + Capture.insert( pair< ZAI, map<ZAI, double> > ( ZAI(95,242,1), toAdd ) ); } map<ZAI, int> index_inver; map<int, ZAI> index; @@ -459,7 +481,7 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is for(it = ZAIDecay.begin() ; it != ZAIDecay.end(); it++) { index.insert( pair<int, ZAI > ( i, (*it).first ) ); - index_inver.insert( pair<ZAI, int > ( (*it).first , i )); + index_inver.insert( pair<ZAI, int > ( (*it).first , i )); i++; } } @@ -468,9 +490,9 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is for(int i = 0; i < (int)index.size(); i++) for(int j = 0; j < (int)index.size(); j++) DecayMatrix[i][j] = 0; - - -// Fill the Decay Part of the Bateman Matrix + + + // Fill the Decay Part of the Bateman Matrix { int i = 0; map<ZAI, pair<double, map< ZAI, double > > >::iterator it; @@ -480,17 +502,17 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is map< ZAI, double > decaylist = (*it).second.second; for(it2 = decaylist.begin(); it2!= decaylist.end(); it2++) { - + map<ZAI, int >::iterator it3 = index_inver.find( (*it2).first ); if( it3 != index_inver.end() ) DecayMatrix[(*it3).second][i] = log(2.)/(*it).second.first * (*it2).second; else { map<ZAI, map<ZAI, double> >::iterator it4 = FastDecay.find( (*it2).first ); - - if( it4 == FastDecay.end() ) + + if( it4 == FastDecay.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it2).first.Z() << " " << (*it2).first.A() << " " << (*it2).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it2).first.Z() << " " << (*it2).first.A() << " " << (*it2).first.I() << endl; exit(1); } @@ -501,25 +523,24 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is it3 = index_inver.find( (*it5).first ); if( it3 == index_inver.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it2).first.Z() << " " << (*it2).first.A() << " " << (*it2).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it2).first.Z() << " " << (*it2).first.A() << " " << (*it2).first.I() << endl; exit(1); } DecayMatrix[(*it3).second][i] = log(2.)/(*it).second.first * (*it2).second * (*it5).second; } - + } } DecayMatrix[i][i] += -log(2.)/(*it).second.first; i++; - + } } - - - double NormFactor = 1; - - + + + + //-------------------------// //--- Perform Evolution ---// //-------------------------// @@ -529,17 +550,17 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is { // Filling the t=0 State; map<ZAI, double > isotopicquantity = isotopicvector.GetIsotopicQuantity(); TMatrixT<double> N_0Matrix = TMatrixT<double>( index.size(),1) ; - + map<ZAI, double >::iterator it ; for(int i = 0; i < (int)index.size(); i++) N_0Matrix[i] = 0; - + for(it = isotopicquantity.begin(); it != isotopicquantity.end(); it++) { map<ZAI, int >::iterator it2; - if( (*it).first.Z() < 90 ) + if( (*it).first.Z() < 90 ) it2 = index_inver.find( ZAI(-2,-2,-2) ); else it2 = index_inver.find( (*it).first ); @@ -547,46 +568,76 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is it2 = index_inver.find( ZAI(-3,-3,-3) ); N_0Matrix[ (*it2).second ][0] = (*it).second ; - + } NMatrix.push_back(N_0Matrix); } - + for(int i = 0; i < 16; i++) { + + double TStep = cycletime/16*i; - + TMatrixT<double> BatemanMatrix = TMatrixT<double>(index.size(),index.size()); BatemanMatrix = DecayMatrix ; - - + + IsotopicVector IVStep; for(int k = 0; k < (int)index.size(); k++) IVStep += index.find(k)->second * NMatrix.back()[k][0]; - + EvolutiveProduct evolutiveproductStep = GetDistancesTo(IVStep, TStep).begin()->second; - double Flux = evolutiveproductStep.GetFlux()->Eval(TStep)*Power/evolutiveproductStep.GetPower(); + double NormFactor = 1; + { + IsotopicVector WantedHMIV = isotopicvector.GetSpeciesComposition(90) + + isotopicvector.GetSpeciesComposition(92) + + isotopicvector.GetSpeciesComposition(93) + + isotopicvector.GetSpeciesComposition(94) + + isotopicvector.GetSpeciesComposition(95) + + isotopicvector.GetSpeciesComposition(96); + + IsotopicVector DBHMIV = evolutiveproductStep.GetIsotopicVectorAt(0).GetSpeciesComposition(90) + + evolutiveproductStep.GetIsotopicVectorAt(0).GetSpeciesComposition(92) + + evolutiveproductStep.GetIsotopicVectorAt(0).GetSpeciesComposition(93) + + evolutiveproductStep.GetIsotopicVectorAt(0).GetSpeciesComposition(94) + + evolutiveproductStep.GetIsotopicVectorAt(0).GetSpeciesComposition(95) + + evolutiveproductStep.GetIsotopicVectorAt(0).GetSpeciesComposition(96); + + NormFactor = Norme(WantedHMIV)/ Norme(DBHMIV); + } + + + double Flux = evolutiveproductStep.GetFlux()->Eval(TStep)*Power/(evolutiveproductStep.GetPower()*NormFactor); + map<ZAI ,TGraph* >::iterator it; -// ---------------- A(n,.) X+Y - for(it = evolutiveproductStep.GetFissionXS().begin() ; it != evolutiveproductStep.GetFissionXS().end(); it++) + // ---------------- A(n,.) X+Y + + map<ZAI ,TGraph* > FissionXS = evolutiveproductStep.GetFissionXS(); + + for(it = FissionXS.begin() ; it != FissionXS.end(); it++) { - if( index_inver.find( (*it).first ) != index_inver.end() ) + + if( index_inver.find( (*it).first ) != index_inver.end() ) { double y; y = (*it).second->Eval(TStep); + BatemanMatrix[ index_inver.find( (*it).first )->second ][index_inver.find( (*it).first )->second] += -y* 1e-24 *Flux; BatemanMatrix[1][ index_inver.find( (*it).first )->second] += 2*y* 1e-24 *Flux; } + } -// ---------------- A(n,.)A+1 - for(it = evolutiveproductStep.GetCaptureXS().begin() ; it != evolutiveproductStep.GetCaptureXS().end(); it++) + // ---------------- A(n,.)A+1 + map<ZAI ,TGraph* > CaptureXS = evolutiveproductStep.GetCaptureXS(); + for(it = CaptureXS.begin(); it != CaptureXS.end(); it++) { - - if( index_inver.find( (*it).first ) != index_inver.end() ) + + if( index_inver.find( (*it).first ) != index_inver.end() ) { double y; y = (*it).second->Eval(TStep); @@ -604,18 +655,18 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is if( it6 != index_inver.end() ) { BatemanMatrix[(*it6).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux ; - + } else { map<ZAI, map<ZAI, double> >::iterator it4 = FastDecay.find( ZAI( (*it).first.Z(), (*it).first.A()+1, (*it).first.I()) ); - - if( it4 == FastDecay.end() ) + + if( it4 == FastDecay.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it).first.Z() << " " << (*it).first.A()+1 << " " << (*it).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it).first.Z() << " " << (*it).first.A()+1 << " " << (*it).first.I() << endl; exit(1); } - + map< ZAI, double >::iterator it5; map< ZAI, double > decaylist2 = (*it4).second; for(it5 = decaylist2.begin(); it5!= decaylist2.end(); it5++) @@ -624,64 +675,65 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is it6 = index_inver.find( (*it5).first ); if( it6 == index_inver.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it).first.Z() << " " << (*it).first.A() << " " << (*it).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it).first.Z() << " " << (*it).first.A() << " " << (*it).first.I() << endl; exit(1); } - + BatemanMatrix[(*it6).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux * (*it5).second; } } } else { - // if( (*it3).first.Z() == 90 && (*it3).first.A() == 232) cout << y* 1e-24 *Flux << endl; + //if( (*it3).first.Z() == 90 && (*it3).first.A() == 232) cout << y* 1e-24 *Flux << endl; map<ZAI, double>::iterator it4; map<ZAI, double> CaptureList = (*it3).second; for(it4 = CaptureList.begin(); it4 != CaptureList.end() ; it4++) { - + map<ZAI, int >::iterator it6 = index_inver.find( (*it4).first ); - - + + if( it6 != index_inver.end() ) BatemanMatrix[(*it6).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux * (*it4).second ; else { map<ZAI, map<ZAI, double> >::iterator it7 = FastDecay.find( (*it4).first ); - - if( it7 == FastDecay.end() ) + + if( it7 == FastDecay.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it7).first.Z() << " " << (*it7).first.A() << " " << (*it7).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it7).first.Z() << " " << (*it7).first.A() << " " << (*it7).first.I() << endl; exit(1); } - + map< ZAI, double >::iterator it5; map< ZAI, double > decaylist2 = (*it7).second; for(it5 = decaylist2.begin(); it5!= decaylist2.end(); it5++) { - + it6 = index_inver.find( (*it5).first ); if( it6 == index_inver.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it7).first.Z() << " " << (*it7).first.A() << " " << (*it7).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it7).first.Z() << " " << (*it7).first.A() << " " << (*it7).first.I() << endl; exit(1); } -// if( (*it6).first.Z() == 92 && (*it6).first.A() == 233) cout << y* 1e-24 *Flux * (*it5).second << endl; + //if( (*it6).first.Z() == 92 && (*it6).first.A() == 233) cout << y* 1e-24 *Flux * (*it5).second << endl; BatemanMatrix[(*it6).second][index_inver.find( (*it).first )->second] += y * 1e-24 * Flux * (*it5).second * (*it4).second; } } - + } } - - + + } } - -// ---------------- A(n,2n)A-1 - for(it = evolutiveproductStep.Getn2nXS().begin() ; it != evolutiveproductStep.Getn2nXS().end(); it++) + + // ---------------- A(n,2n)A-1 + map<ZAI ,TGraph* > n2nXS = evolutiveproductStep.Getn2nXS(); + for(it = n2nXS.begin() ; it != n2nXS.end(); it++) { - if( index_inver.find( (*it).first ) != index_inver.end() ) + if( index_inver.find( (*it).first ) != index_inver.end() ) { double y; y = (*it).second->Eval(TStep); @@ -690,14 +742,14 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is map<ZAI, int>::iterator it3 = index_inver.find( ZAI( (*it).first.Z(), (*it).first.A()-1, 0) ); - if( it3 != index_inver.end() ) - BatemanMatrix[(*it3).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux; + if( it3 != index_inver.end() ) + BatemanMatrix[(*it3).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux; else { - + map<ZAI, map<ZAI, double> >::iterator it4 = FastDecay.find( ZAI( (*it).first.Z(), (*it).first.A()-1, 0) ); - - if( it4 == FastDecay.end() ) + + if( it4 == FastDecay.end() ) { it3 = index_inver.find( ZAI( -3, -3, -3 ) ); BatemanMatrix[(*it3).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux; @@ -708,11 +760,11 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is map< ZAI, double > decaylist2 = (*it4).second; for(it5 = decaylist2.begin(); it5!= decaylist2.end(); it5++) { - + it3 = index_inver.find( (*it5).first ); if( it3 == index_inver.end() ) { - cout << "Problem in FastDecay for nuclei " << (*it4).first.Z() << " " << (*it4).first.A() << " " << (*it4).first.I() << endl; + cout << "Problem in FastDecay for nuclei " << (*it4).first.Z() << " " << (*it4).first.A() << " " << (*it4).first.I() << endl; exit(1); } BatemanMatrix[(*it3).second][index_inver.find( (*it).first )->second] += y* 1e-24 *Flux * (*it5).second ; @@ -721,10 +773,10 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is } } } - + // ---------------- Evolution - TMatrixT<double> NEvoluingMatrix = TMatrixT<double>(index.size(),1); - + TMatrixT<double> NEvolutionMatrix = TMatrixT<double>(index.size(),1); + double TStepMax = cycletime/16.; timevector[i+1] = timevector[i] + TStepMax; @@ -737,40 +789,40 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is if(k == j) IdMatrix[j][k] = 1; else IdMatrix[j][k] = 0; } - - + + TMatrixT<double> BatemanMatrixDL = TMatrixT<double>(index.size(),index.size()); // Order 0 Term from the DL : Id TMatrixT<double> BatemanMatrixDLTermN = TMatrixT<double>(index.size(),index.size()); // Addind it; - - + + { BatemanMatrixDLTermN = IdMatrix; BatemanMatrixDL = BatemanMatrixDLTermN; - - + + int j = 1; double NormN = 0; - + do { NormN = 0; - TMatrixT<double> BatemanMatrixDLTermtmp = TMatrixT<double>(index.size(),index.size()); // Addind it; + TMatrixT<double> BatemanMatrixDLTermtmp = TMatrixT<double>(index.size(),index.size()); // Adding it; BatemanMatrixDLTermtmp = BatemanMatrixDLTermN; BatemanMatrixDLTermN.Mult(BatemanMatrixDLTermtmp, BatemanMatrix ); - + BatemanMatrixDLTermN *= 1./j; BatemanMatrixDL += BatemanMatrixDLTermN; - + NormN = 0; for(int m = 0; m < (int)index.size(); m++) for(int n = 0; n < (int)index.size(); n++) - NormN += BatemanMatrixDLTermN[m][n]*BatemanMatrixDLTermN[m][n]; + NormN += BatemanMatrixDLTermN[m][n]*BatemanMatrixDLTermN[m][n]; j++; } while ( NormN != 0); } - NEvoluingMatrix = BatemanMatrixDL * NMatrix.back() ; - NMatrix.push_back(NEvoluingMatrix); + NEvolutionMatrix = BatemanMatrixDL * NMatrix.back() ; + NMatrix.push_back(NEvolutionMatrix); } @@ -781,15 +833,15 @@ EvolutiveProduct EvolutionDataBase<IsotopicVector>::GenerateDB(IsotopicVector is double ZAIQuantity[NMatrix.size()]; for(int j = 0; j < (int)NMatrix.size(); j++) ZAIQuantity[j] = (NMatrix[j])[i][0]; - + GeneratedDB.Insert(pair<ZAI, TGraph*> (index.find(i)->second, new TGraph(NMatrix.size(), timevector, ZAIQuantity) ) ); } - GeneratedDB.SetPower(Power*NormFactor ); + GeneratedDB.SetPower(Power ); GeneratedDB.SetFuelType(fFuelType ); -// GeneratedDB.SetReactorType(fReactorType ); -// GeneratedDB.SetHMMass(fHMMass*NormFactor ); - + //GeneratedDB.SetReactorType(fReactorType ); + //GeneratedDB.SetHMMass(fHMMass*NormFactor ); + return GeneratedDB; DBGL; } diff --git a/source/trunk/src/EvolutiveProduct.cxx b/source/trunk/src/EvolutiveProduct.cxx index f0558ed79..be1e2df13 100755 --- a/source/trunk/src/EvolutiveProduct.cxx +++ b/source/trunk/src/EvolutiveProduct.cxx @@ -408,7 +408,7 @@ void EvolutiveProduct::ReadDB(string DBfile) start = 0; getline(InfoDB, line); if (StringLine::NextWord(line, start, ' ') == "CycleTime") - fReactorType = StringLine::NextWord(line, start, ' '); + fCycleTime = atof(StringLine::NextWord(line, start, ' ').c_str());; getline(InfoDB, line); // Assembly HM Mass start = 0; getline(InfoDB, line); @@ -948,7 +948,7 @@ EvolutiveProduct EvolutiveProduct::GenerateDBFor(IsotopicVector isotopicvector) } // ---------------- Evolution - TMatrixT<double> NEvoluingMatrix = TMatrixT<double>(index.size(),1); + TMatrixT<double> NEvolutionMatrix = TMatrixT<double>(index.size(),1); double TStepMax; { @@ -985,7 +985,7 @@ EvolutiveProduct EvolutiveProduct::GenerateDBFor(IsotopicVector isotopicvector) { NormN = 0; - TMatrixT<double> BatemanMatrixDLTermtmp = TMatrixT<double>(index.size(),index.size()); // Addind it; + TMatrixT<double> BatemanMatrixDLTermtmp = TMatrixT<double>(index.size(),index.size()); // Adding it; BatemanMatrixDLTermtmp = BatemanMatrixDLTermN; BatemanMatrixDLTermN.Mult(BatemanMatrixDLTermtmp, BatemanMatrix ); @@ -1002,8 +1002,8 @@ EvolutiveProduct EvolutiveProduct::GenerateDBFor(IsotopicVector isotopicvector) } - NEvoluingMatrix = BatemanMatrixDL * NMatrix.back() ; - NMatrix.push_back(NEvoluingMatrix); + NEvolutionMatrix = BatemanMatrixDL * NMatrix.back() ; + NMatrix.push_back(NEvolutionMatrix); } diff --git a/source/trunk/src/FabricationPlant.cxx b/source/trunk/src/FabricationPlant.cxx index 2243bc941..9e67c81f5 100644 --- a/source/trunk/src/FabricationPlant.cxx +++ b/source/trunk/src/FabricationPlant.cxx @@ -21,14 +21,14 @@ #include <cmath> #include <algorithm> -//________________________________________________________________________ -// -// FabricationPlant -// -// -// -// -//________________________________________________________________________ + //________________________________________________________________________ + // + // FabricationPlant + // + // + // + // + //________________________________________________________________________ ClassImp(FabricationPlant) template <class T> T random(T a, T b) //peak random numebr between a and b @@ -41,35 +41,60 @@ template <class T> T random(T a, T b) //peak random numebr between a and b FabricationPlant::FabricationPlant() { -DBGL; + DBGL; fChronologicalTimePriority = false; -DBGL; + fFabricationTime = -1; + fUpdateReferenceDBatEachStep = true; + + + // Warning + + cout << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + cout << "\t Chronological Stock Priority set! "<< endl << endl; + + fLog->fLog << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + fLog->fLog <<"\t Chronological Stock Priority set! "<< endl << endl; + + DBGL; } FabricationPlant::FabricationPlant(Storage* storage, Storage* reusable, double fabircationtime) { -DBGL; - fFabricationTime = (cSecond)fabircationtime; + DBGL; + fChronologicalTimePriority = false; + fUpdateReferenceDBatEachStep = true; + + + fFabricationTime = (cSecond)fabircationtime; fStorage = storage; fReUsable = reusable; -DBGL; + + + cout << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + cout << "\t Chronological Stock Priority has been set! "<< endl; + cout << "\t Fabrication time set to \t " << (double)(fFabricationTime/3600/24/365.25) << " year" << endl << endl; + + fLog->fLog << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + fLog->fLog <<"\t Chronological Stock Priority has been set! "<< endl; + fLog->fLog << "\t Fabrication time set to \t " << (double)(fFabricationTime/3600/24/365.25) << " year" << endl << endl; + DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ FabricationPlant::~FabricationPlant() { -DBGL; -DBGL; + DBGL; + DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ void FabricationPlant::AddValorisableIV(ZAI zai, double factor) { -DBGL; + DBGL; pair<map<ZAI, double>::iterator, bool> IResult; if(factor > 1) factor = 1; @@ -79,14 +104,14 @@ DBGL; if(IResult.second == false) IResult.first->second = factor; } -DBGL; + DBGL; } -//________________________________________________________________________ -//_______________________________ Evolution ______________________________ -//________________________________________________________________________ + //________________________________________________________________________ + //_______________________________ Evolution ______________________________ + //________________________________________________________________________ void FabricationPlant::Evolution(cSecond t) { DBGL; @@ -99,7 +124,7 @@ void FabricationPlant::Evolution(cSecond t) DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ void FabricationPlant::FabricationPlantEvolution(cSecond t) { DBGL; @@ -207,11 +232,11 @@ void FabricationPlant::BuildFuelForReactor(int ReactorId) { map<ZAI ,double>::iterator it; map<ZAI ,double> isotopicquantity = GetDecay( stock , fFabricationTime).GetSpeciesComposition(94).GetIsotopicQuantity(); - + for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++ ) { if ((*it).first.A() >= 238 && (*it).first.A() <= 242) - { + { nPu_1 += (*it).second; Sum_AlphaI_nPuI += FuelType->GetFuelParameter()[(*it).first.A() -237]*(*it).second; } @@ -219,21 +244,21 @@ void FabricationPlant::BuildFuelForReactor(int ReactorId) isotopicquantity = stock.GetSpeciesComposition(94).GetIsotopicQuantity(); for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++ ) - if ((*it).first.A() >= 238 && (*it).first.A() <= 242) - { - MPu_1 += (*it).second * (ZAImass.find( (*it).first )->second)/Na*1e-6; - } + if ((*it).first.A() >= 238 && (*it).first.A() <= 242) + { + MPu_1 += (*it).second * (ZAImass.find( (*it).first )->second)/Na*1e-6; + } isotopicquantity = GetDecay( FullUsedStock , fFabricationTime).GetSpeciesComposition(94).GetIsotopicQuantity(); for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++ ) if ((*it).first.A() >= 238 && (*it).first.A() <= 242) - { + { Sum_AlphaI_nPuI0 += FuelType->GetFuelParameter()[(*it).first.A() -237]*(*it).second; - } + } } double StockFactionToUse = 0; - + double NT = HMmass*1e6 * Na / ZAImass.find( ZAI(92,238,0) )->second; double N1 = (BU - FuelType->GetFuelParameter()[6]) * NT; @@ -307,9 +332,9 @@ void FabricationPlant::BuildFuelForReactor(int ReactorId) -//________________________________________________________________________ -//_________________________________ Decay ________________________________ -//________________________________________________________________________ + //________________________________________________________________________ + //_________________________________ Decay ________________________________ + //________________________________________________________________________ IsotopicVector FabricationPlant::GetDecay(IsotopicVector isotopicvector, cSecond t) { DBGL; @@ -331,26 +356,35 @@ IsotopicVector FabricationPlant::GetDecay(IsotopicVector isotopicvector, cSecond } -//________________________________________________________________________ -//_____________________________ Reactor & DB _____________________________ -//________________________________________________________________________ + //________________________________________________________________________ + //_____________________________ Reactor & DB _____________________________ + //________________________________________________________________________ EvolutiveProduct FabricationPlant::BuildEvolutiveDB(int ReactorId,IsotopicVector isotopicvector) { -DBGL; + DBGL; EvolutionDataBase<IsotopicVector>* evolutiondb = fParc->GetReactor()[ReactorId]->GetFuelType(); - - isotopicvector = GetDecay(isotopicvector, fFabricationTime); - map<double, EvolutiveProduct> distances = evolutiondb->GetDistancesTo(isotopicvector); + isotopicvector = GetDecay(isotopicvector, fFabricationTime); - - EvolutiveProduct EvolBuild = distances.begin()->second.GenerateDBFor(isotopicvector); + EvolutiveProduct EvolBuild; + if( fUpdateReferenceDBatEachStep == true ) + { + EvolutiveProduct EvolBuild = evolutiondb->GenerateDB(isotopicvector, + fParc->GetReactor()[ReactorId]->GetCycleTime(), + fParc->GetReactor()[ReactorId]->GetPower()); + } + else + { + map<double, EvolutiveProduct> distances = evolutiondb->GetDistancesTo(isotopicvector); + EvolBuild = distances.begin()->second.GenerateDBFor(isotopicvector); + + } return EvolBuild; -DBGL; + DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ void FabricationPlant::TakeReactorFuel(int Id) { DBGL; @@ -362,7 +396,7 @@ void FabricationPlant::TakeReactorFuel(int Id) DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ EvolutiveProduct FabricationPlant::GetReactorEvolutionDB(int ReactorId) { DBGL; @@ -371,9 +405,9 @@ EvolutiveProduct FabricationPlant::GetReactorEvolutionDB(int ReactorId) DBGL; } -//________________________________________________________________________ -//_______________________________ Storage ________________________________ -//________________________________________________________________________ + //________________________________________________________________________ + //_______________________________ Storage ________________________________ + //________________________________________________________________________ IsotopicVector FabricationPlant::GetStockToRecycle() { @@ -396,30 +430,30 @@ IsotopicVector FabricationPlant::GetStockToRecycle() DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ void FabricationPlant::RecycleStock(double fraction) { -DBGL; + DBGL; fFractionToTake.back().second = fraction; -DBGL; + DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ void FabricationPlant::DumpStock() { -DBGL; + DBGL; -DBGL; + DBGL; } -//________________________________________________________________________ + //________________________________________________________________________ pair<IsotopicVector, IsotopicVector> FabricationPlant::Separation(IsotopicVector isotopicvector) { -DBGL; - //[0] = re-use ; [1] = waste + DBGL; + //[0] = re-use ; [1] = waste pair<IsotopicVector, IsotopicVector> IVTmp; map<ZAI ,double> isotopicquantity = isotopicvector.GetIsotopicQuantity(); @@ -435,7 +469,7 @@ DBGL; } else IVTmp.second.Add( (*it).first, (*it).second ); //waste } -DBGL; + DBGL; return IVTmp; } diff --git a/source/trunk/src/Reactor.cxx b/source/trunk/src/Reactor.cxx index 8397b61f3..501e78b79 100755 --- a/source/trunk/src/Reactor.cxx +++ b/source/trunk/src/Reactor.cxx @@ -56,6 +56,21 @@ Reactor::Reactor(EvolutionDataBase<IsotopicVector>* fueltypeDB, fCreationTime = (cSecond)creationtime; fLifeTime = (cSecond)lifetime; + + cout << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + cout << "\t Fuel Composition is not fixed ! "<< endl; + cout << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + cout << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl << endl; + cout << "!!WARNING!! !!!Reactor!!! You need to set Burn-up/Power/CycleTime (2 of 3) & Heavy Metal Mass Manualy !! " << endl; + + + fLog->fLog << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + fLog->fLog << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + fLog->fLog << "\t Fuel Composition is not fixed ! "<< endl; + fLog->fLog << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + fLog->fLog << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl << endl; DBGL; } @@ -87,6 +102,31 @@ Reactor::Reactor(double Power, EvolutionDataBase<IsotopicVector>* fueltypeDB, fCreationTime = (cSecond)creationtime; fLifeTime = (cSecond)lifetime; + + + + cout << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + cout << "\t Fuel Composition is not fixed ! "<< endl; + cout << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + cout << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t The Effective Thermal Power set at \t " << (double)(fPower *1e-6) << " MW" << endl; + cout << "\t Burn-Up at end of Cycle set at \t " << (double)(fBurnUp) << " GWj/t" << endl; + cout << "\t The corresponding Cycle Time is\t " << (double)(fCycleTime/3600/24/365.25) << " year" << endl; + cout << "\t The Heavy Metal Mass in the Core set at " << (double)(fHeavyMetalMass) << " tons" << endl << endl; + + + fLog->fLog << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + fLog->fLog << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + fLog->fLog << "\t Fuel Composition is not fixed ! "<< endl; + fLog->fLog << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + fLog->fLog << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t The Effective Thermal Power set at \t " << (double)(fPower *1e-6) << " MW" << endl; + fLog->fLog << "\t Burn-Up at end of Cycle set at \t " << (double)(fBurnUp) << " GWj/t" << endl; + fLog->fLog << "\t The corresponding Cycle Time is\t " << (double)(fCycleTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t The Heavy Metal Mass in the Core set at " << (double)(fHeavyMetalMass) << " tons" << endl << endl; + DBGL; } @@ -117,6 +157,30 @@ Reactor::Reactor(EvolutionDataBase<IsotopicVector>* fueltypeDB, fCreationTime = (cSecond)creationtime; fLifeTime = (cSecond)lifetime; fPower = BurnUp*3600.*24. / (fCycleTime) * HMMass *1e9; //BU in GWd/t + + + cout << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + cout << "\t Fuel Composition is not fixed ! "<< endl; + cout << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + cout << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t The Cycle Time set at\t " << (double)(fCycleTime/3600/24/365.25) << " year" << endl; + cout << "\t Burn-Up at end of Cycle set at \t " << (double)(fBurnUp) << " GWj/t" << endl; + cout << "\t The corresponding Effective Thermal Power is \t " << (double)(fPower *1e-6) << " MW" << endl; + cout << "\t The Heavy Metal Mass in the Core set at " << (double)(fHeavyMetalMass) << " tons" << endl << endl; + + + fLog->fLog << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + fLog->fLog << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + fLog->fLog << "\t Fuel Composition is not fixed ! "<< endl; + fLog->fLog << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + fLog->fLog << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t The Cycle Time set at\t " << (double)(fCycleTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t Burn-Up at end of Cycle set at \t " << (double)(fBurnUp) << " GWj/t" << endl; + fLog->fLog << "\t The corresponding Effective Thermal Power is \t " << (double)(fPower *1e-6) << " MW" << endl; + fLog->fLog << "\t The Heavy Metal Mass in the Core set at " << (double)(fHeavyMetalMass) << " tons" << endl << endl; + DBGL; } @@ -149,6 +213,28 @@ Reactor::Reactor(EvolutiveProduct evolutivedb, fIVBeginCycle = fEvolutionDB.GetIsotopicVectorAt(0); fIVInCycle = fEvolutionDB.GetIsotopicVectorAt(0); fIVOutCycle = fEvolutionDB.GetIsotopicVectorAt( (cSecond)(fCycleTime/fEvolutionDB.GetPower()*fPower) ); + + + cout << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + cout << "\t Fuel Composition is fixed ! "<< endl; + cout << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + cout << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + cout << "\t The Cycle Time set at\t " << (double)(fCycleTime/3600/24/365.25) << " year" << endl; + cout << "\t The Effective Thermal Power is \t " << (double)(fPower *1e-6) << " MW" << endl; + cout << "\t The Heavy Metal Mass in the Core set at " << (double)(fHeavyMetalMass) << " tons" << endl << endl; + + + fLog->fLog << "!!Info!! !!!FabricationPlant!!! A FabricationPlant has been define :" << endl; + fLog->fLog << "!!Info!! !!!Reactor!!! A Reactor has been define :" << endl; + fLog->fLog << "\t Fuel Composition is not fixed ! "<< endl; + fLog->fLog << "\t Fuel Type set to : \t "<< fFuelTypeDB->GetFuelType() << endl; + fLog->fLog << "\t Creation time set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t Life time (Operating's Duration) set at \t " << (double)(fCreationTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t The Cycle Time set at\t " << (double)(fCycleTime/3600/24/365.25) << " year" << endl; + fLog->fLog << "\t The Effective Thermal Power is \t " << (double)(fPower *1e-6) << " MW" << endl; + fLog->fLog << "\t The Heavy Metal Mass in the Core set at " << (double)(fHeavyMetalMass) << " tons" << endl << endl; + DBGL; } -- GitLab