Skip to content
Snippets Groups Projects
Commit a214edf9 authored by BaM's avatar BaM
Browse files

updating source: step 3 : removing the CLASSV3 Branches 

git-svn-id: svn+ssh://svn.in2p3.fr/class@479 0e7d625b-0364-4367-a6be-d5be4a48d228
parent 05226e12
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 0 additions and 3174 deletions
source/branches/CLASSV3/Model/Equivalence/EQM_LIN_PWR_MOX.cxx deleted 100644 → 0
+ 0
− 262
View file @ 05226e12
#include "EQM_LIN_PWR_MOX.hxx"
#include "CLASSConstante.hxx"
#include <vector>
#include "StringLine.hxx"
#include "CLASSLogger.hxx"
#include "IsotopicVector.hxx"
EQM_LIN_PWR_MOX::EQM_LIN_PWR_MOX(string WeightPath):EquivalenceModel(new CLASSLogger("EQM_LIN_PWR_MOX.log"))
{
fWeightPath = WeightPath;
ifstream DataDB(fWeightPath.c_str()); // Open the File
if(!DataDB)
WARNING << "Can't open \"" << fWeightPath << "\"\n" << endl;
string line;
int start = 0; // First Get Fuel Parameter
getline(DataDB, line);
if( StringLine::NextWord(line, start, ' ') != "PARAM")
{
ERROR << " Bad Database file : " << fWeightPath << " Can't find the Parameter of the DataBase " << endl;
exit (1);
}
while(start < (int)line.size())
fFuelParameter.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
INFO << " " << (int)fFuelParameter.size() << " parameters have been read " << endl;
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
// ADD ENrichment of the U reading !!!!!!!!!!!!!!!!!!!!!!!!!!!! //
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
ZAI U8(92,238,0);
ZAI U5(92,235,0);
double U5_enrich= 0.0025;
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
ZAI Pu8(94,238,0);
ZAI Pu9(94,239,0);
ZAI Pu0(94,240,0);
ZAI Pu1(94,241,0);
ZAI Pu2(94,242,0);
fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
}
EQM_LIN_PWR_MOX::EQM_LIN_PWR_MOX(CLASSLogger* log, string WeightPath):EquivalenceModel(log)
{
fWeightPath = WeightPath;
ifstream DataDB(fWeightPath.c_str()); // Open the File
if(!DataDB)
WARNING << " Can't open \"" << fWeightPath << "\"\n" << endl;
string line;
int start = 0; // First Get Fuel Parameter
getline(DataDB, line);
if( StringLine::NextWord(line, start, ' ') != "PARAM")
{
ERROR << " Bad Database file : " << fWeightPath << " Can't find the Parameter of the DataBase"<< endl;
exit (1);
}
while(start < (int)line.size())
fFuelParameter.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
INFO << fFuelParameter.size() << " have been read"<< endl;
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
// ADD ENrichment of the U reading !!!!!!!!!!!!!!!!!!!!!!!!!!!! //
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
ZAI U8(92,238,0);
ZAI U5(92,235,0);
double U5_enrich= 0.0025;
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
ZAI Pu8(94,238,0);
ZAI Pu9(94,239,0);
ZAI Pu0(94,240,0);
ZAI Pu1(94,241,0);
ZAI Pu2(94,242,0);
fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
}
EQM_LIN_PWR_MOX::~EQM_LIN_PWR_MOX()
{
}
//________________________________________________________________________
vector<double> EQM_LIN_PWR_MOX::BuildFuel(double BurnUp, double HMMass,vector<IsotopicVector> FissilArray, vector<IsotopicVector> FertilArray)
{
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
// ADD ENrichment of the U check ++ Check Un seul fertile !!!! //
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
//-----------------------------------------------------------------------------//
vector<double> lambda;
for(int i = 0; i < (int) (FissilArray.size() + FertilArray.size()); i++)
lambda.push_back(0);
double Na = 6.02214129e23; //N Avogadro
IsotopicVector FullUsedStock;
IsotopicVector stock;
bool FuelBuild = false;
if(FissilArray.size() == 0)
{
for(int i = 0; i < (int)lambda.size(); i++)
lambda[i] = -1;
FuelBuild = true;
}
int N_FissilStock_OnCheck = 0;
while(!FuelBuild)
{
double nPu_0 = 0;
double MPu_0 = 0;
{
map<ZAI ,double>::iterator it;
map<ZAI ,double> isotopicquantity = FullUsedStock.GetSpeciesComposition(94).GetIsotopicQuantity();
for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++ )
nPu_0 += (*it).second;
IsotopicVector IV_Pm_Am = (FullUsedStock.GetSpeciesComposition(94)
+ ZAI(94,241,0)*FullUsedStock.GetZAIIsotopicQuantity(95,241,0));
//Add the 241Am as 241Pu... the Pu is not old in the Eq Model but is in the FissileArray....;
MPu_0 += cZAIMass.GetMass(IV_Pm_Am);
}
stock = FissilArray[N_FissilStock_OnCheck];
double nPu_1 = 0;
double MPu_1 = 0;
double Sum_AlphaI_nPuI = 0;
double Sum_AlphaI_nPuI0 = 0;
{
map<ZAI ,double>::iterator it;
map<ZAI ,double> isotopicquantity = stock.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 += fFuelParameter[(*it).first.A() -237]*(*it).second;
}
}
isotopicquantity = (stock.GetSpeciesComposition(94) + ZAI(94,241,0)*stock.GetZAIIsotopicQuantity(95,241,0)).GetIsotopicQuantity();
IsotopicVector IV_Pm_Am = (stock.GetSpeciesComposition(94)
+ ZAI(94,241,0)*stock.GetZAIIsotopicQuantity(95,241,0));
//Add the 241Am as 241Pu... the Pu is not old in the Eq Model but is in the FissileArray....
MPu_1 += cZAIMass.GetMass(IV_Pm_Am);
isotopicquantity = FullUsedStock.GetSpeciesComposition(94).GetIsotopicQuantity();
for( it = isotopicquantity.begin(); it != isotopicquantity.end(); it++ )
if ((*it).first.A() >= 238 && (*it).first.A() <= 242)
{
Sum_AlphaI_nPuI0 += fFuelParameter[(*it).first.A() -237]*(*it).second;
}
}
double StockFactionToUse = 0;
double NT = HMMass*1e6 * Na / (cZAIMass.GetMass( ZAI(92,238,0) ) * 0.997
+ cZAIMass.GetMass( ZAI(92,235,0) ) * 0.003 );
double N1 = (BurnUp - fFuelParameter[6]) * NT;
double N2 = -Sum_AlphaI_nPuI0;
double N3 = -fFuelParameter[0] * Na / (cZAIMass.GetMass( ZAI(92,238,0) )*0.997
+ cZAIMass.GetMass( ZAI(92,235,0) )*0.003 )
* (HMMass*1e6 - MPu_0*1e6);
double D1 = Sum_AlphaI_nPuI;
double D2 = -fFuelParameter[0] * MPu_1*1e6 * Na / (cZAIMass.GetMass( ZAI(92,238,0) )*0.997
+ cZAIMass.GetMass( ZAI(92,235,0) )*0.003 );
StockFactionToUse = (N1 + N2 + N3) / (D1 + D2);
if(StockFactionToUse < 0)
{
WARNING << "!!!FabricationPlant!!! Oups Bug in calculating stock fraction to use "<< endl;
lambda[N_FissilStock_OnCheck] = 0.;
N_FissilStock_OnCheck++;
FuelBuild = false;
}
else if( StockFactionToUse > 1 )
{
FullUsedStock += stock;
lambda[N_FissilStock_OnCheck] = 1;
N_FissilStock_OnCheck++;
FuelBuild = false;
}
else
{
lambda[N_FissilStock_OnCheck] = StockFactionToUse;
FuelBuild = true;
double U8_Quantity = (HMMass - (MPu_0+StockFactionToUse*MPu_1 ))/(cZAIMass.GetMass( ZAI(92,238,0))*0.997 + cZAIMass.GetMass( ZAI(92,235,0))*0.003 )*Na/1e-6;
lambda.back() = U8_Quantity / FertilArray[0].GetSumOfAll();
}
if( N_FissilStock_OnCheck == (int) FissilArray.size() ) // Check if the last Fissil stock has been tested... quit if so...
{
for(int i = 0; i < (int)lambda.size(); i++)
lambda[i] = -1;
FuelBuild = true;
}
}
return lambda;
}
source/branches/CLASSV3/Model/Equivalence/EQM_LIN_PWR_MOX.hxx deleted 100644 → 0
+ 0
− 48
View file @ 05226e12
#ifndef _EQM_LIN_PWR_MOX_HXX
#define _EQM_LIN_PWR_MOX_HXX
#include "EquivalenceModel.hxx"
#include <string>
/*!
\file
\brief Header file for EQM_MLP_MOX class.
@author BaM
@version 1.0
*/
using namespace std;
//-----------------------------------------------------------------------------//
/*!
Define a EQM_MLP_MOX.
The aim of these class is to constuct a fuel from an equivalence model
based on a Linear Eq Model BU = SUM (alpha_i*n_i}
@author BaM
@version 3.0
*/
//________________________________________________________________________
class EQM_LIN_PWR_MOX : public EquivalenceModel
{
public :
EQM_LIN_PWR_MOX(string WeightPath);
EQM_LIN_PWR_MOX(CLASSLogger* log, string WeightPath);
~EQM_LIN_PWR_MOX();
virtual vector<double> BuildFuel(double BurnUp, double HMMass, vector<IsotopicVector> FissilArray, vector<IsotopicVector> FertilArray );
private :
string fWeightPath;
vector<double> fFuelParameter;
};
#endif
source/branches/CLASSV3/Model/Equivalence/EQM_MLP_PWR_MOX.cxx deleted 100644 → 0
+ 0
− 182
View file @ 05226e12
#include "EquivalenceModel.hxx"
#include "EQM_MLP_PWR_MOX.hxx"
#include "CLASSLogger.hxx"
#include "StringLine.hxx"
#include <string>
#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
#include <cassert>
#include "TSystem.h"
#include "TMVA/Reader.h"
#include "TMVA/Tools.h"
#include "TMVA/MethodCuts.h"
//________________________________________________________________________
//
// EQM_MLP_MOX
//
// Equivalenve Model based on multi layer perceptron from TMVA (root cern)
// For REP MOX use
//
//________________________________________________________________________
EQM_MLP_MOX::EQM_MLP_MOX(string TMVAWeightPath):EquivalenceModel(new CLASSLogger("EQM_MLP_MOX.log"))
{
fTMVAWeightPath = TMVAWeightPath;
ZAI U8(92,238,0);
ZAI U5(92,235,0);
double U5_enrich= 0.0025;
ZAI Pu8(94,238,0);
ZAI Pu9(94,239,0);
ZAI Pu0(94,240,0);
ZAI Pu1(94,241,0);
ZAI Pu2(94,242,0);
fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
INFO<<"__An equivalence model of PWR MOX has been define__"<<endl;
INFO<<"\tThis model is based on a multi layer perceptron"<<endl;
INFO<<"\t\tThe TMVA weight file is :"<<endl;
INFO<<"\t\t\t"<<fTMVAWeightPath<<endl;
}
//________________________________________________________________________
EQM_MLP_MOX::EQM_MLP_MOX(CLASSLogger* log, string TMVAWeightPath):EquivalenceModel(log)
{
fTMVAWeightPath = TMVAWeightPath;
ZAI U8(92,238,0);
ZAI U5(92,235,0);
double U5_enrich= 0.0025;
ZAI Pu8(94,238,0);
ZAI Pu9(94,239,0);
ZAI Pu0(94,240,0);
ZAI Pu1(94,241,0);
ZAI Pu2(94,242,0);
fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
INFO<<"__An equivalence model of PWR MOX has been define__"<<endl;
INFO<<"\tThis model is based on a multi layer perceptron"<<endl;
INFO<<"\t\tThe TMVA weight file is :"<<endl;
INFO<<"\t\t\t"<<fTMVAWeightPath<<endl;
}
//________________________________________________________________________
TTree* EQM_MLP_MOX::CreateTMVAInputTree(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp)
{
TTree* InputTree = new TTree("EQTMP", "EQTMP");
float Pu8 = 0;
float Pu9 = 0;
float Pu10 = 0;
float Pu11 = 0;
float Pu12 = 0;
float Am1 = 0;
float U5_enrichment = 0;
float BU = 0;
InputTree->Branch( "Pu8" ,&Pu8 ,"Pu8/F" );
InputTree->Branch( "Pu9" ,&Pu9 ,"Pu9/F" );
InputTree->Branch( "Pu10" ,&Pu10 ,"Pu10/F" );
InputTree->Branch( "Pu11" ,&Pu11 ,"Pu11/F" );
InputTree->Branch( "Pu12" ,&Pu12 ,"Pu12/F" );
InputTree->Branch( "Am1" ,&Am1 ,"Am1/F" );
InputTree->Branch( "U5_enrichment" ,&U5_enrichment ,"U5_enrichment/F" );
InputTree->Branch( "BU" ,&BU ,"BU/F" );
float U8 = Fertil.GetZAIIsotopicQuantity(92,238,0);
float U5 = Fertil.GetZAIIsotopicQuantity(92,235,0);
float U4 = Fertil.GetZAIIsotopicQuantity(92,234,0);
float UTOT = U8 + U5 + U4;
Pu8 = Fissil.GetZAIIsotopicQuantity(94,238,0);
Pu9 = Fissil.GetZAIIsotopicQuantity(94,239,0);
Pu10 = Fissil.GetZAIIsotopicQuantity(94,240,0);
Pu11 = Fissil.GetZAIIsotopicQuantity(94,241,0);
Pu12 = Fissil.GetZAIIsotopicQuantity(94,242,0);
Am1 = Fissil.GetZAIIsotopicQuantity(95,241,0);
double TOTPU=(Pu8+Pu9+Pu10+Pu11+Pu12+Am1);
Pu8 = Pu8 / TOTPU;
Pu9 = Pu9 / TOTPU;
Pu10= Pu10 / TOTPU;
Pu11= Pu11 / TOTPU;
Pu12= Pu12 / TOTPU;
Am1 = Am1 / TOTPU;
U5_enrichment = U5 / UTOT;
BU=BurnUp;
//cout<<"Pu8 "<<Pu8 <<" Pu9 "<< Pu9 <<" Pu10 "<< Pu10 << " Pu11 "<< Pu11 <<" Pu12 "<<Pu12 <<" Am1 "<<Am1<<endl;
//cout<<"BU "<<BU<<" U5_enrichment "<<U5_enrichment<<endl;
if(Pu8 + Pu9 + Pu10 + Pu11 + Pu12 + Am1 > 1.00001 )//?????1.00001??? I don't know it! goes in condition if =1 !! may be float/double issue ...
{
ERROR << Pu8 << " " << Pu9 << " " << Pu10 << " " << Pu11 << " " << Pu12 << " " << Am1 << endl;
exit(0);
}
// All value are molar (!weight)
InputTree->Fill();
return InputTree;
}
//________________________________________________________________________
double EQM_MLP_MOX::ExecuteTMVA(TTree* theTree)
{
// --- Create the Reader object
TMVA::Reader *reader = new TMVA::Reader( "Silent" );
// Create a set of variables and declare them to the reader
// - the variable names MUST corresponds in name and type to those given in the weight file(s) used
Float_t Pu8,Pu9,Pu10,Pu11,Pu12,Am1,BU,U5_enrichment;
reader->AddVariable( "BU" ,&BU );
reader->AddVariable( "U5_enrichment",&U5_enrichment );
reader->AddVariable( "Pu8" ,&Pu8 );
reader->AddVariable( "Pu9" ,&Pu9 );
reader->AddVariable( "Pu10" ,&Pu10);
reader->AddVariable( "Pu11" ,&Pu11);
reader->AddVariable( "Pu12" ,&Pu12);
reader->AddVariable( "Am1" ,&Am1 );
// --- Book the MVA methods
// Book method MLP
TString methodName = "MLP method";
reader->BookMVA( methodName, fTMVAWeightPath );
//theTree->SetBranchAddress("teneur",&teneur);
theTree->SetBranchAddress( "BU" ,&BU );
theTree->SetBranchAddress( "U5_enrichment" ,&U5_enrichment ) ;
theTree->SetBranchAddress( "Pu8" ,&Pu8 );
theTree->SetBranchAddress( "Pu9" ,&Pu9 );
theTree->SetBranchAddress( "Pu10" ,&Pu10 );
theTree->SetBranchAddress( "Pu11" ,&Pu11 );
theTree->SetBranchAddress( "Pu12" ,&Pu12 );
theTree->SetBranchAddress( "Am1" ,&Am1 );
theTree->GetEntry(0);
Float_t val = (reader->EvaluateRegression( methodName ))[0];
delete reader;
delete theTree;
return (double)val; //retourne teneur
}
//________________________________________________________________________
double EQM_MLP_MOX::GetFissileMolarFraction(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp)
{DBGL
return ExecuteTMVA(CreateTMVAInputTree(Fissil,Fertil,BurnUp));
}
source/branches/CLASSV3/Model/Equivalence/EQM_MLP_PWR_MOX.hxx deleted 100644 → 0
+ 0
− 50
View file @ 05226e12
#ifndef _EQM_MLP_MOX_HXX
#define _EQM_MLP_MOX_HXX
#include "EquivalenceModel.hxx"
#include "TTree.h"
/*!
\file
\brief Header file for EQM_MLP_MOX class.
@author BaM
@version 1.0
*/
using namespace std;
//-----------------------------------------------------------------------------//
/*!
Define a EQM_MLP_MOX.
The aim of these class is to constuct a fuel from an equivalence model
based on a Multi layer perceptron
@author BaM
@version 3.0
*/
//________________________________________________________________________
class EQM_MLP_MOX : public EquivalenceModel
{
public :
EQM_MLP_MOX(string TMVAWeightPath); //!< Constructor string TMVAWeightPath => PATH/TMVAWeight.xml (path to tmva weight)
EQM_MLP_MOX(CLASSLogger* log, string TMVAWeightPath); //!< Constructor CLASSLogger* log ,string TMVAWeightPath => PATH/TMVAWeight.xml
virtual double GetFissileMolarFraction(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp); //!<Return the molar fraction of fissile element thanks to a Multi Layer Perceptron
private :
TTree* CreateTMVAInputTree(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp);//!<Create input tmva tree to be read by ExecuteTMVA
double ExecuteTMVA(TTree* theTree);//!<Execute the MLP according to the input tree created by CreateTMVAInputTree
string fTMVAWeightPath;;//!<The weight needed by TMVA to construct and execute the multilayer perceptron
};
#endif
source/branches/CLASSV3/Model/Equivalence/EQM_QUAD_PWR_MOX.cxx deleted 100644 → 0
+ 0
− 136
View file @ 05226e12
#include "EQM_QUAD_PWR_MOX.hxx"
#include <vector>
#include "StringLine.hxx"
#include "CLASSLogger.hxx"
EQM_QUAD_PWR_MOX::EQM_QUAD_PWR_MOX(string WeightPath):EquivalenceModel(new CLASSLogger("EQM_QUAD_PWR_MOX.log"))
{
fWeightPath = WeightPath;
ifstream DataDB(fWeightPath.c_str()); // Open the File
if(!DataDB)
WARNING << " Can't open \"" << fWeightPath << "\"" << endl;
string line;
int start = 0; // First Get Fuel Parameter
getline(DataDB, line);
if( StringLine::NextWord(line, start, ' ') != "PARAM")
{
ERROR << "Bad Database file : " << fWeightPath << " Can't find the Parameter of the DataBase " << endl;
exit (1);
}
while(start < (int)line.size())
fFuelParameter.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
INFO << fFuelParameter.size() << " have been read " << endl;
ZAI U8(92,238,0);
ZAI U5(92,235,0);
double U5_enrich= 0.0025;
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
ZAI Pu8(94,238,0);
ZAI Pu9(94,239,0);
ZAI Pu0(94,240,0);
ZAI Pu1(94,241,0);
ZAI Pu2(94,242,0);
fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
}
EQM_QUAD_PWR_MOX::EQM_QUAD_PWR_MOX(CLASSLogger* log, string WeightPath):EquivalenceModel(log)
{
fWeightPath = WeightPath;
ifstream DataDB(fWeightPath.c_str()); // Open the File
if(!DataDB)
WARNING << " Can't open \"" << fWeightPath << "\"" << endl;
string line;
int start = 0; // First Get Fuel Parameter
getline(DataDB, line);
if( StringLine::NextWord(line, start, ' ') != "PARAM")
{
ERROR << "Bad Database file : " << fWeightPath << " Can't find the Parameter of the DataBase " << endl;
exit (1);
}
while(start < (int)line.size())
fFuelParameter.push_back(atof(StringLine::NextWord(line, start, ' ').c_str()));
INFO << fFuelParameter.size() << " have been read " << endl;
ZAI U8(92,238,0);
ZAI U5(92,235,0);
double U5_enrich= 0.0025;
fFertileList = U5*U5_enrich + U8*(1-U5_enrich);
ZAI Pu8(94,238,0);
ZAI Pu9(94,239,0);
ZAI Pu0(94,240,0);
ZAI Pu1(94,241,0);
ZAI Pu2(94,242,0);
fFissileList = Pu8*1+Pu9*1+Pu0*1+Pu1*1+Pu2*1;
}
EQM_QUAD_PWR_MOX::~EQM_QUAD_PWR_MOX()
{
}
double EQM_QUAD_PWR_MOX::GetFissileMolarFraction(IsotopicVector Fissile,IsotopicVector Fertile,double BurnUp)
{
ZAI ZAIList[6] = {ZAI(94,238,0), ZAI(94,239,0), ZAI(94,240,0), ZAI(94,241,0), ZAI(94,242,0), ZAI(95,241,0) };
vector<double> PuCompo;
double Sum = Fissile.GetSumOfAll();
for(int i = 0; i< 5; i++)
PuCompo.push_back( Fissile.GetZAIIsotopicQuantity(ZAIList[i])/Sum);
PuCompo[2] += Fissile.GetZAIIsotopicQuantity(ZAIList[5])/Sum;
double A = 0;
if(PuCompo[0] <= PuCompo[2] && PuCompo[0] <= PuCompo[4] && PuCompo[1] + PuCompo[3] >= 0.40 && PuCompo[1] >0 )
{
int par = 0;
for(int j = 0 ; j < 5 ; j++)
{
A += fFuelParameter[par] * PuCompo[j] ;
par++;
for(int i = j ; i < 5 ; i++)
{
A += fFuelParameter[par] *PuCompo[i] *PuCompo[j];
par++;
}
}
A += fFuelParameter[par];
}
else
{
cout << "the composition is not in the range of the Model" << endl;
}
return A;
}
source/branches/CLASSV3/Model/Equivalence/EQM_QUAD_PWR_MOX.hxx deleted 100644 → 0
+ 0
− 48
View file @ 05226e12
#ifndef _EQM_QUAD_PWR_MOX_HXX
#define _EQM_QUAD_PWR_MOX_HXX
#include "EquivalenceModel.hxx"
#include <string>
/*!
\file
\brief Header file for EQM_MLP_MOX class.
@author BaM
@version 1.0
*/
using namespace std;
//-----------------------------------------------------------------------------//
/*!
Define a EQM_MLP_MOX.
The aim of these class is to constuct a fuel from an equivalence model
based on a Quadratic Pu equivalent Model
@author BaM
@version 3.0
*/
//________________________________________________________________________
class EQM_QUAD_PWR_MOX : public EquivalenceModel
{
public :
EQM_QUAD_PWR_MOX(string WeightPath);
EQM_QUAD_PWR_MOX(CLASSLogger* log, string WeightPath);
~EQM_QUAD_PWR_MOX();
virtual double GetFissileMolarFraction(IsotopicVector Fissil,IsotopicVector Fertil,double BurnUp);
private :
string fWeightPath;
vector<double> fFuelParameter;
};
#endif
source/branches/CLASSV3/Model/Irradiation/IM_Matrix.cxx deleted 100644 → 0
+ 0
− 301
View file @ 05226e12
//
// IM_Matrix.cxx
// CLASSSource
//
// Created by BaM on 04/05/2014.
// Copyright (c) 2014 BaM. All rights reserved.
//
#include "IM_Matrix.hxx"
#include "IsotopicVector.hxx"
#include "CLASSConstante.hxx"
#include "CLASSLogger.hxx"
#include "StringLine.hxx"
#include <TGraph.h>
#include <TString.h>
#include <sstream>
#include <string>
#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
using namespace std;
//________________________________________________________________________
IM_Matrix::IM_Matrix():IrradiationModel(new CLASSLogger("IM_Matrix.log"))
{
fShorstestHalflife = 3600.*24*160.; //cut by default all nuclei with a shorter liftime than the Cm242 -> remain 33 actinides
}
IM_Matrix::IM_Matrix(CLASSLogger* log):IrradiationModel(log)
{
fShorstestHalflife = 3600.*24*160.; //cut by default all nuclei with a shorter liftime than the Cm242 -> remain 33 actinides
}
//________________________________________________________________________
/* Evolution Calculation */
//________________________________________________________________________
EvolutionData IM_Matrix::GenerateEvolutionData(IsotopicVector isotopicvector, EvolutionData XSSet, double Power, double cycletime)
{
DBGL
if(fFastDecay.size() == 0)
NuclearDataInitialization();
string ReactorType;
vector< TMatrixT<double> > NMatrix ;// TMatrixT<double>(decayindex.size(),1))
{ // Filling the t=0 State;
map<ZAI, double > isotopicquantity = isotopicvector.GetIsotopicQuantity();
TMatrixT<double> N_0Matrix = TMatrixT<double>( fReverseMatrixIndex.size(),1) ;
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
N_0Matrix[i] = 0;
map<ZAI, double >::iterator it ;
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
N_0Matrix[i] = 0;
for(it = isotopicquantity.begin(); it != isotopicquantity.end(); it++)
{
/// Need TO change with FP managment
map<ZAI, int >::iterator it2;
if( (*it).first.Z() < fZAIThreshold )
it2 = fMatrixIndex.find( ZAI(-2,-2,-2) );
else it2 = fMatrixIndex.find( (*it).first );
if(it2 == fMatrixIndex.end() ) //If not in index should be TMP, can't be fast decay for new Fuel !!!
it2 = fMatrixIndex.find( ZAI(-3,-3,-3) );
N_0Matrix[ (*it2).second ][0] = (*it).second ;
}
isotopicquantity.clear();
NMatrix.push_back(N_0Matrix);
N_0Matrix.Clear();
}
//-------------------------//
//--- Perform Evolution ---//
//-------------------------//
ReactorType = XSSet.GetReactorType();
double M_ref = XSSet.GetHeavyMetalMass();
double M = cZAIMass.GetMass(isotopicvector.GetActinidesComposition());
double Power_ref = XSSet.GetPower();
int NStep = XSSet.GetFissionXS().begin()->second->GetN();
double* DBTimeStep = XSSet.GetFissionXS().begin()->second->GetX();
int InsideStep = 10;
double timevector[NStep];
timevector[0] = 0;
double Flux[NStep];
TMatrixT<double> FissionEnergy = TMatrixT<double>(fReverseMatrixIndex.size(),1);
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
FissionEnergy[i] = 0;
{
map< ZAI, int >::iterator it;
for(it = fMatrixIndex.begin(); it != fMatrixIndex.end(); it++)
{
map< ZAI, double >::iterator it2 = fFissionEnergy.find(it->first);
if(it2 == fFissionEnergy.end())
{
if(it->first.Z() > fZAIThreshold)
FissionEnergy[it->second][0] = 1.9679e6*it->first.A()-2.601e8; // //simple linear fit to known values ;extrapolation to unknown isotopes
else FissionEnergy[it->second][0] = 0;
}
else
FissionEnergy[it->second][0] = it2->second;
}
}
vector< TMatrixT<double> > FissionXSMatrix; // Store The Fisison XS Matrix
vector< TMatrixT<double> > CaptureXSMatrix; // Store The Capture XS Matrix
vector< TMatrixT<double> > n2nXSMatrix; // Store The n2N XS Matrix
DBGL
for(int i = 0; i < NStep-1; i++)
{
double TStepMax = ( (DBTimeStep[i+1]-DBTimeStep[i] ) ) * Power_ref/M_ref / Power*M ; // Get the next Time step
TMatrixT<double> BatemanMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size());
TMatrixT<double> BatemanReactionMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size());
TMatrixT<double> NEvolutionMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),1);
NEvolutionMatrix = NMatrix.back();
FissionXSMatrix.push_back(GetFissionXsMatrix(XSSet, DBTimeStep[i])); //Feel the fission reaction Matrix
CaptureXSMatrix.push_back(GetCaptureXsMatrix(XSSet, DBTimeStep[i])); //Feel the capture reaction Matrix
n2nXSMatrix.push_back(Getn2nXsMatrix(XSSet, DBTimeStep[i])); //Feel the (n,2n) reaction Matrix
// ---------------- Evolution
BatemanReactionMatrix = FissionXSMatrix[i];
BatemanReactionMatrix += CaptureXSMatrix[i];
BatemanReactionMatrix += n2nXSMatrix[i];
for(int k=0; k < InsideStep; k++)
{
double ESigmaN = 0;
for (int j = 0; j < (int)fReverseMatrixIndex.size() ; j++)
ESigmaN -= FissionXSMatrix[i][j][j]*NEvolutionMatrix[j][0]*1.6e-19*FissionEnergy[j][0];
// Update Flux
double Flux_k = Power/ESigmaN;
if(k==0)
Flux[i]=Flux_k;
BatemanMatrix = BatemanReactionMatrix;
BatemanMatrix *= Flux_k;
BatemanMatrix += fDecayMatrix ;
BatemanMatrix *= TStepMax/InsideStep ;
TMatrixT<double> IdMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size());
for(int j = 0; j < (int)fReverseMatrixIndex.size(); j++)
for(int k = 0; k < (int)fReverseMatrixIndex.size(); k++)
{
if(k == j) IdMatrix[j][k] = 1;
else IdMatrix[j][k] = 0;
}
TMatrixT<double> BatemanMatrixDL = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size()); // Order 0 Term from the DL : Id
TMatrixT<double> BatemanMatrixDLTermN = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size()); // Addind it;
{
BatemanMatrix *= TStepMax ;
BatemanMatrixDLTermN = IdMatrix;
BatemanMatrixDL = BatemanMatrixDLTermN;
int j = 1;
double NormN;
do
{
TMatrixT<double> BatemanMatrixDLTermtmp = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size()); // Adding it;
BatemanMatrixDLTermtmp = BatemanMatrixDLTermN;
BatemanMatrixDLTermN.Mult(BatemanMatrixDLTermtmp, BatemanMatrix );
BatemanMatrixDLTermN *= 1./j;
BatemanMatrixDL += BatemanMatrixDLTermN;
NormN = 0;
for(int m = 0; m < (int)fReverseMatrixIndex.size(); m++)
for(int n = 0; n < (int)fReverseMatrixIndex.size(); n++)
NormN += BatemanMatrixDLTermN[m][n]*BatemanMatrixDLTermN[m][n];
j++;
} while ( NormN != 0 );
}
NEvolutionMatrix = BatemanMatrixDL * NEvolutionMatrix ;
}
NMatrix.push_back(NEvolutionMatrix);
timevector[i+1] = timevector[i] + TStepMax;
BatemanMatrix.Clear();
BatemanReactionMatrix.Clear();
NEvolutionMatrix.Clear();
}
DBGL
FissionXSMatrix.push_back(GetFissionXsMatrix(XSSet, DBTimeStep[NStep-1])); //Feel the reaction Matrix
CaptureXSMatrix.push_back(GetCaptureXsMatrix(XSSet, DBTimeStep[NStep-1])); //Feel the reaction Matrix
n2nXSMatrix.push_back(Getn2nXsMatrix(XSSet, DBTimeStep[NStep-1])); //Feel the reaction Matrix
EvolutionData GeneratedDB = EvolutionData(GetLog());
double ESigmaN = 0;
for (int j = 0; j < (int)fReverseMatrixIndex.size() ; j++)
ESigmaN -= FissionXSMatrix.back()[j][j]*NMatrix.back()[j][0]*1.6e-19*FissionEnergy[j][0];
Flux[NStep-1] = Power/ESigmaN;
GeneratedDB.SetFlux( new TGraph(NStep, timevector, Flux) );
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
{
double ZAIQuantity[NMatrix.size()];
double FissionXS[NStep];
double CaptureXS[NStep];
double n2nXS[NStep];
for(int j = 0; j < (int)NMatrix.size(); j++)
ZAIQuantity[j] = (NMatrix[j])[i][0];
for(int j = 0; j < NStep; j++)
{
FissionXS[j] = FissionXSMatrix[j][i][i];
CaptureXS[j] = CaptureXSMatrix[j][i][i];
n2nXS[j] = n2nXSMatrix[j][i][i];
}
GeneratedDB.NucleiInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NMatrix.size(), timevector, ZAIQuantity)));
GeneratedDB.FissionXSInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NStep, timevector, FissionXS)));
GeneratedDB.CaptureXSInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NStep, timevector, CaptureXS)));
GeneratedDB.n2nXSInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NStep, timevector, n2nXS)));
}
GeneratedDB.SetPower(Power );
GeneratedDB.SetHeavyMetalMass(M);
GeneratedDB.SetReactorType(ReactorType );
GeneratedDB.SetCycleTime(cycletime);
for (int i = 0; i < (int) FissionXSMatrix.size(); i++)
{
FissionXSMatrix[i].Clear();
CaptureXSMatrix[i].Clear();
n2nXSMatrix[i].Clear();
}
FissionXSMatrix.clear();
CaptureXSMatrix.clear();
n2nXSMatrix.clear();
DBGL
return GeneratedDB;
}
source/branches/CLASSV3/Model/Irradiation/IM_Matrix.hxx deleted 100644 → 0
+ 0
− 63
View file @ 05226e12
#ifndef _IMMATRIX_HXX
#define _IMMATRIX_HXX
/*!
\file
\brief Header file for IrradiationModel class.
@author BaM
@version 2.0
*/
#include "IrradiationModel.hxx"
using namespace std;
class CLASSLogger;
//-----------------------------------------------------------------------------//
/*!
Define a IM_Matrix.
The aim of these class is perform the calculation of the evolution of a fuel trough irradiation solving numericaly the Bateman using Matrix.
@author BaM
@version 3.0
*/
//________________________________________________________________________
class EvolutionData;
class IM_Matrix : public IrradiationModel
{
public :
IM_Matrix();
IM_Matrix(CLASSLogger* log);
/// virtueal method called to perform the irradiation calculation using a set of cross section.
/*!
Perform the Irradiation Calcultion using the XSSet data
\param IsotopicVector IV isotopic vector to irradiate
\param EvolutionData XSSet set of corss section to use to perform the evolution calculation
*/
virtual EvolutionData GenerateEvolutionData(IsotopicVector IV, EvolutionData XSSet, double Power, double cycletime);
//}
private :
};
#endif
source/branches/CLASSV3/Model/Irradiation/IM_RK4.cxx deleted 100644 → 0
+ 0
− 408
View file @ 05226e12
//
// IM_RK4.cxx
// CLASSSource
//
// Created by BaM on 04/05/2014.
// Copyright (c) 2014 BaM. All rights reserved.
//
#include "IM_RK4.hxx"
#include "IsotopicVector.hxx"
#include "CLASSConstante.hxx"
#include "CLASSLogger.hxx"
#include "StringLine.hxx"
#include <TGraph.h>
#include <TString.h>
#include <sstream>
#include <string>
#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
using namespace std;
//________________________________________________________________________
IM_RK4::IM_RK4():IrradiationModel(new CLASSLogger("IM_RK4.log")), DynamicalSystem()
{
fTheNucleiVector = 0;
fTheMatrix = 0;
SetForbidNegativeValue();
}
IM_RK4::IM_RK4(CLASSLogger* log):IrradiationModel(log), DynamicalSystem()
{
fTheNucleiVector = 0;
fTheMatrix = 0;
SetForbidNegativeValue();
}
//________________________________________________________________________
/* Evolution Calculation */
//________________________________________________________________________
EvolutionData IM_RK4::GenerateEvolutionData(IsotopicVector isotopicvector, EvolutionData XSSet, double Power, double cycletime)
{
DBGL
if(fFastDecay.size() == 0)
{
NuclearDataInitialization();
fNVar = fReverseMatrixIndex.size();
}
SetTheMatrixToZero();
SetTheNucleiVectorToZero();
string ReactorType;
vector< TMatrixT<double> > NMatrix ;// TMatrixT<double>(decayindex.size(),1))
{ // Filling the t=0 State;
map<ZAI, double > isotopicquantity = isotopicvector.GetIsotopicQuantity();
TMatrixT<double> N_0Matrix = TMatrixT<double>( fReverseMatrixIndex.size(),1) ;
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
N_0Matrix[i] = 0;
map<ZAI, double >::iterator it ;
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
N_0Matrix[i] = 0;
for(it = isotopicquantity.begin(); it != isotopicquantity.end(); it++)
{
/// Need TO change with FP managment
map<ZAI, int >::iterator it2;
if( (*it).first.Z() < fZAIThreshold )
it2 = fMatrixIndex.find( ZAI(-2,-2,-2) );
else it2 = fMatrixIndex.find( (*it).first );
if(it2 == fMatrixIndex.end() ) //If not in index should be TMP, can't be fast decay for new Fuel !!!
it2 = fMatrixIndex.find( ZAI(-3,-3,-3) );
N_0Matrix[ (*it2).second ][0] = (*it).second ;
}
isotopicquantity.clear();
NMatrix.push_back(N_0Matrix);
N_0Matrix.Clear();
}
//-------------------------//
//--- Perform Evolution ---//
//-------------------------//
ReactorType = XSSet.GetReactorType();
double M_ref = XSSet.GetHeavyMetalMass();
double M = cZAIMass.GetMass(isotopicvector.GetActinidesComposition());
double Power_ref = XSSet.GetPower();
int NStep = XSSet.GetFissionXS().begin()->second->GetN();
double* DBTimeStep = XSSet.GetFissionXS().begin()->second->GetX();
int InsideStep = 10;
double timevector[NStep];
timevector[0] = 0;
double Flux[NStep];
TMatrixT<double> FissionEnergy = TMatrixT<double>(fReverseMatrixIndex.size(),1);
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
FissionEnergy[i] = 0;
{
map< ZAI, int >::iterator it;
for(it = fMatrixIndex.begin(); it != fMatrixIndex.end(); it++)
{
map< ZAI, double >::iterator it2 = fFissionEnergy.find(it->first);
if(it2 == fFissionEnergy.end())
{
if(it->first.Z() > fZAIThreshold)
FissionEnergy[it->second][0] = 1.9679e6*it->first.A()-2.601e8; // //simple linear fit to known values ;extrapolation to unknown isotopes
else FissionEnergy[it->second][0] = 0;
}
else
FissionEnergy[it->second][0] = it2->second;
}
}
vector< TMatrixT<double> > FissionXSMatrix; // Store The Fisison XS Matrix
vector< TMatrixT<double> > CaptureXSMatrix; // Store The Capture XS Matrix
vector< TMatrixT<double> > n2nXSMatrix; // Store The n2N XS Matrix
DBGL
for(int i = 0; i < NStep-1; i++)
{
double TStepMax = ( (DBTimeStep[i+1]-DBTimeStep[i] ) ) * Power_ref/M_ref / Power*M ; // Get the next Time step
TMatrixT<double> BatemanMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size());
TMatrixT<double> BatemanReactionMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size());
TMatrixT<double> NEvolutionMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),1);
NEvolutionMatrix = NMatrix.back();
FissionXSMatrix.push_back(GetFissionXsMatrix(XSSet, DBTimeStep[i])); //Feel the fission reaction Matrix
CaptureXSMatrix.push_back(GetCaptureXsMatrix(XSSet, DBTimeStep[i])); //Feel the capture reaction Matrix
n2nXSMatrix.push_back(Getn2nXsMatrix(XSSet, DBTimeStep[i])); //Feel the (n,2n) reaction Matrix
// ---------------- Evolution
BatemanReactionMatrix = FissionXSMatrix[i];
/* for(int u=0;u<fReverseMatrixIndex.size();u++)
for(int v=0;v<fReverseMatrixIndex.size();v++)
if(CaptureXSMatrix[i][u][v]!=0)
cout<<CaptureXSMatrix[i][u][v]<<endl;
exit(0);
*/
BatemanReactionMatrix += CaptureXSMatrix[i];
BatemanReactionMatrix += n2nXSMatrix[i];
for(int k=0; k < InsideStep; k++)
{
double ESigmaN = 0;
for (int j = 0; j < (int)fReverseMatrixIndex.size() ; j++)
ESigmaN -= FissionXSMatrix[i][j][j]*NEvolutionMatrix[j][0]*1.6e-19*FissionEnergy[j][0];
// Update Flux
double Flux_k = Power/ESigmaN;
if(k==0)
Flux[i]=Flux_k;
BatemanMatrix = BatemanReactionMatrix;
BatemanMatrix *= Flux_k;
BatemanMatrix += fDecayMatrix ;
/*for(int u=0;u<fReverseMatrixIndex.size();u++)
for(int v=0;v<fReverseMatrixIndex.size();v++)
if(BatemanMatrix[u][v]!=0)
cout<<BatemanMatrix[u][v]<<endl;
exit(1);*/
SetTheMatrixToZero();
SetTheNucleiVectorToZero();
SetTheMatrix(BatemanMatrix);
SetTheNucleiVector(NEvolutionMatrix);
RungeKutta(fTheNucleiVector, timevector[i]+TStepMax/InsideStep*k, timevector[i]+TStepMax/InsideStep*(k+1), fNVar);
NEvolutionMatrix = GetTheNucleiVector();
}
NEvolutionMatrix = GetTheNucleiVector();
NMatrix.push_back(NEvolutionMatrix);
timevector[i+1] = timevector[i] + TStepMax;
BatemanMatrix.Clear();
BatemanReactionMatrix.Clear();
NEvolutionMatrix.Clear();
}
FissionXSMatrix.push_back(GetFissionXsMatrix(XSSet, DBTimeStep[NStep-1])); //Feel the reaction Matrix
CaptureXSMatrix.push_back(GetCaptureXsMatrix(XSSet, DBTimeStep[NStep-1])); //Feel the reaction Matrix
n2nXSMatrix.push_back(Getn2nXsMatrix(XSSet, DBTimeStep[NStep-1])); //Feel the reaction Matrix
EvolutionData GeneratedDB = EvolutionData(GetLog());
double ESigmaN = 0;
for (int j = 0; j < (int)fReverseMatrixIndex.size() ; j++)
ESigmaN -= FissionXSMatrix.back()[j][j]*NMatrix.back()[j][0]*1.6e-19*FissionEnergy[j][0];
Flux[NStep-1] = Power/ESigmaN;
GeneratedDB.SetFlux( new TGraph(NStep, timevector, Flux) );
for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++)
{
double ZAIQuantity[NMatrix.size()];
double FissionXS[NStep];
double CaptureXS[NStep];
double n2nXS[NStep];
for(int j = 0; j < (int)NMatrix.size(); j++)
ZAIQuantity[j] = (NMatrix[j])[i][0];
for(int j = 0; j < NStep; j++)
{
FissionXS[j] = FissionXSMatrix[j][i][i];
CaptureXS[j] = CaptureXSMatrix[j][i][i];
n2nXS[j] = n2nXSMatrix[j][i][i];
}
GeneratedDB.NucleiInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NMatrix.size(), timevector, ZAIQuantity)));
GeneratedDB.FissionXSInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NStep, timevector, FissionXS)));
GeneratedDB.CaptureXSInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NStep, timevector, CaptureXS)));
GeneratedDB.n2nXSInsert(pair<ZAI, TGraph*> (fReverseMatrixIndex[i], new TGraph(NStep, timevector, n2nXS)));
}
DBGL
GeneratedDB.SetPower(Power );
GeneratedDB.SetHeavyMetalMass(M);
GeneratedDB.SetReactorType(ReactorType );
GeneratedDB.SetCycleTime(cycletime);
ResetTheMatrix();
ResetTheNucleiVector();
for (int i = 0; i < (int) FissionXSMatrix.size(); i++)
{
FissionXSMatrix[i].Clear();
CaptureXSMatrix[i].Clear();
n2nXSMatrix[i].Clear();
}
FissionXSMatrix.clear();
CaptureXSMatrix.clear();
n2nXSMatrix.clear();
DBGL
return GeneratedDB;
}
void IM_RK4::ResetTheMatrix()
{
if(fTheMatrix)
{
for(int i= 0; i<fNVar; i++)
delete [] fTheMatrix[i];
delete [] fTheMatrix;
}
fTheMatrix = 0;
}
void IM_RK4::SetTheMatrixToZero()
{
ResetTheMatrix();
fNVar = fReverseMatrixIndex.size();
fTheMatrix = new double*[fNVar];
#pragma omp parallel for
for(int i= 0; i < fNVar; i++)
fTheMatrix[i] = new double[fNVar];
for(int i = 0; i < fNVar; i++)
for(int k = 0; k < fNVar; k++)
{
fTheMatrix[i][k]=0.0;
}
}
//________________________________________________________________________
void IM_RK4::ResetTheNucleiVector()
{
if(fTheNucleiVector)
delete [] fTheNucleiVector;
fTheNucleiVector = 0;
}
//________________________________________________________________________
void IM_RK4::SetTheNucleiVectorToZero()
{
ResetTheNucleiVector();
fTheNucleiVector = new double[fNVar];
#pragma omp parallel for
for(int i = 0; i < fNVar; i++)
fTheNucleiVector[i]=0.0;
}
//________________________________________________________________________
void IM_RK4::BuildEqns(double t, double *N, double *dNdt)
{
double sum=0;
// pragma omp parallel for reduction(+:sum)
for(int i = 0; i < fNVar; i++)
{
sum=0;
for(int k = 0; k < fNVar; k++)
{
sum += fTheMatrix[i][k]*N[k];
}
dNdt[i] = sum;
}
}
//________________________________________________________________________
void IM_RK4::SetTheMatrix(TMatrixT<double> BatemanMatrix)
{
for (int k = 0; k < (int)fNVar; k++)
for (int l = 0; l < (int)fMatrixIndex.size(); l++)
fTheMatrix[l][k] = BatemanMatrix[l][k];
}
//________________________________________________________________________
TMatrixT<double> IM_RK4::GetTheMatrix()
{
TMatrixT<double> BatemanMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),fReverseMatrixIndex.size());
for (int k = 0; k < (int)fNVar; k++)
for (int l = 0; l < (int)fMatrixIndex.size(); l++)
BatemanMatrix[l][k] = fTheMatrix[l][k];
return BatemanMatrix;
}
//________________________________________________________________________
void IM_RK4::SetTheNucleiVector(TMatrixT<double> NEvolutionMatrix)
{
for (int k = 0; k < (int)fNVar; k++)
fTheNucleiVector[k] = NEvolutionMatrix[k][0];
}
//________________________________________________________________________
TMatrixT<double> IM_RK4::GetTheNucleiVector()
{
TMatrixT<double> NEvolutionMatrix = TMatrixT<double>(fReverseMatrixIndex.size(),1);
for (int k = 0; k < (int)fNVar; k++)
NEvolutionMatrix[k][0] = fTheNucleiVector[k];
return NEvolutionMatrix;
}
source/branches/CLASSV3/Model/Irradiation/IM_RK4.hxx deleted 100644 → 0
+ 0
− 102
View file @ 05226e12
#ifndef _IMRK4_HXX
#define _IMRK4_HXX
/*!
\file
\brief Header file for IrradiationModel class.
@author BaM
@version 2.0
*/
#include "DynamicalSystem.hxx"
#include "IrradiationModel.hxx"
using namespace std;
class CLASSLogger;
//-----------------------------------------------------------------------------//
/*!
Define a IM_RK4.
The aim of these class is perform the calculation of the evolution of a fuel trough irradiation solving numericaly the Bateman using RK4.
@author BaM
@version 3.0
*/
//________________________________________________________________________
class EvolutionData;
class IM_RK4 : public IrradiationModel, DynamicalSystem
{
public :
IM_RK4();
IM_RK4(CLASSLogger* Log);
/// virtueal method called to perform the irradiation calculation using a set of cross section.
/*!
Perform the Irradiation Calcultion using the XSSet data
\param IsotopicVector IV isotopic vector to irradiate
\param EvolutionData XSSet set of corss section to use to perform the evolution calculation
*/
virtual EvolutionData GenerateEvolutionData(IsotopicVector IV, EvolutionData XSSet, double Power, double cycletime);
//}
//********* RK4 Method *********//
//@}
/*!
\name RK4 Method
*/
//@{
using DynamicalSystem::RungeKutta;
//! Pre-treatment Runge-Kutta method.
/*!
// This method does initialisation and then call DynamicalSystem::RungeKutta
// \param t1: initial time
// \param t2: final time
*/
void BuildEqns(double t, double *N, double *dNdt);
void SetTheMatrixToZero(); //!< Initialize the evolution Matrix
void ResetTheMatrix();
void SetTheMatrix(TMatrixT<double> BatemanMatrix); //!< Set the Evolution Matrix (Bateman equations)
TMatrixT<double> GetTheMatrix(); //!< return the Evolution Matrix (Bateman equations)
void SetTheNucleiVectorToZero(); //!< Initialize the evolution Matrix
void ResetTheNucleiVector();
void SetTheNucleiVector(TMatrixT<double> NEvolutionMatrix); //!< Set the Evolution Matrix (Bateman equations)
TMatrixT<double> GetTheNucleiVector(); //!< return the Evolution Matrix (Bateman equations)
//@}
private :
double *fTheNucleiVector; //!< The evolving atoms copied from Material proportions.
double **fTheMatrix; //!< The evolution Matrix
double fPrecision; //!< Precision of the RungeKutta
double fHestimate; //!< RK Step estimation.
double fHmin; //!< RK minimum Step.
double fMaxHdid; //!< store the effective RK max step
double fMinHdid; //!< store the effective RK min step
bool fIsNegativeValueAllowed; //!< whether or not negative value are physical.
};
#endif
source/branches/CLASSV3/Model/XS/XSM_CLOSEST.cxx deleted 100644 → 0
+ 0
− 372
View file @ 05226e12
#include "XSModel.hxx"
#include "XSM_CLOSEST.hxx"
#include "CLASSLogger.hxx"
#include "StringLine.hxx"
#include <TGraph.h>
#include <TString.h>
#include "TSystem.h"
#include "TROOT.h"
#include <sstream>
#include <string>
#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
//________________________________________________________________________
//
// XSM_CLOSEST
//
//
//
//
//________________________________________________________________________
XSM_CLOSEST::XSM_CLOSEST(string DB_index_file, bool oldreadmethod ): XSModel(new CLASSLogger("XSM_CLOSEST.log"))
{
fOldReadMethod = oldreadmethod;
fDataBaseIndex = DB_index_file;
fDistanceType = 0;
fWeightedDistance = false;
fEvolutionDataInterpolation = false;
ReadDataBase();
// Warning
INFO << " A EvolutionData has been define :" << endl;
INFO << "\t His index is : \"" << DB_index_file << "\" " << endl;
INFO << "\t " << fFuelDataBank.size() << " EvolutionData have been read."<< endl << endl;
}
//________________________________________________________________________
XSM_CLOSEST::XSM_CLOSEST(CLASSLogger* Log,string DB_index_file, bool oldreadmethod ): XSModel(Log)
{
fOldReadMethod = oldreadmethod;
fDataBaseIndex = DB_index_file;
fDistanceType = 0;
fWeightedDistance = false;
fEvolutionDataInterpolation = false;
ReadDataBase();
// Warning
INFO << " A EvolutionData has been define :" << endl;
INFO << "\t His index is : \"" << DB_index_file << "\" " << endl;
INFO << "\t " << fFuelDataBank.size() << " EvolutionData have been read."<< endl << endl;
}
//________________________________________________________________________
XSM_CLOSEST::~XSM_CLOSEST()
{
for( int i = 0; i < (int)fFuelDataBank.size(); i++)
fFuelDataBank[i].DeleteEvolutionData();
fFuelDataBank.clear();
}
//________________________________________________________________________
void XSM_CLOSEST::ReadDataBase()
{
if(fFuelDataBank.size() != 0)
{
for( int i = 0; i < (int)fFuelDataBank.size(); i++)
fFuelDataBank[i].DeleteEvolutionData();
fFuelDataBank.clear();
}
ifstream DataDB(fDataBaseIndex.c_str()); // Open the File
if(!DataDB)
WARNING << " Can't open \"" << fDataBaseIndex << "\"\n" << endl;
vector<double> vTime;
vector<double> vTimeErr;
string line;
int start = 0;
// First Get Fuel Type
getline(DataDB, line);
if( StringLine::NextWord(line, start, ' ') != "TYPE")
{
ERROR << " Bad Database file : " << fDataBaseIndex << " Can't find the type of the DataBase"<< endl;
exit (1);
}
fFuelType = StringLine::NextWord(line, start, ' ');
//Then Get All the Database
while (!DataDB.eof())
{
getline(DataDB, line);
if(line != "")
{
EvolutionData evolutionproduct(GetLog(), line, fOldReadMethod);
fFuelDataBank.push_back(evolutionproduct);
}
}
}
//________________________________________________________________________
//________________________________________________________________________
//________________________________________________________________________
/* Distance Calculation */
//________________________________________________________________________
//________________________________________________________________________
//________________________________________________________________________
map<double, int> XSM_CLOSEST::GetDistancesTo(IsotopicVector isotopicvector, double t)
{
map<double, int> distances;
for( int i = 0; i < (int)fFuelDataBank.size(); i++)
{
pair<map<double, int>::iterator, bool> IResult;
IsotopicVector ActinidesCompositionAtT = fFuelDataBank[i].GetIsotopicVectorAt(t).GetActinidesComposition();
IsotopicVector IV_ActinidesComposition = isotopicvector.GetActinidesComposition();
double NormalisationFactor = Norme(IV_ActinidesComposition) / Norme( ActinidesCompositionAtT );
double distance = Distance( IV_ActinidesComposition,
ActinidesCompositionAtT / NormalisationFactor,
fDistanceType,
fDistanceParameter);
IResult = distances.insert( pair< double, int >(distance, i) );
}
return distances;
}
//________________________________________________________________________
EvolutionData XSM_CLOSEST::GetCrossSections(IsotopicVector isotopicvector, double t)
{
DBGL
double distance = 0;
int N_BestEvolutionData = 0;
if(fWeightedDistance)
{
DBGL
IsotopicVector ActinidesCompositionAtT = fFuelDataBank[0].GetIsotopicVectorAt(t).GetActinidesComposition();
IsotopicVector IV_ActinidesComposition = isotopicvector.GetActinidesComposition();
double NormalisationFactor = Norme( ActinidesCompositionAtT ) / Norme(IV_ActinidesComposition);
distance = Distance( IV_ActinidesComposition / NormalisationFactor,
fFuelDataBank[0]);
for( int i = 1; i < (int)fFuelDataBank.size(); i++)
{
double D = 0;
ActinidesCompositionAtT = fFuelDataBank[i].GetIsotopicVectorAt(t).GetActinidesComposition();
IV_ActinidesComposition = isotopicvector.GetActinidesComposition();
NormalisationFactor = Norme( ActinidesCompositionAtT ) / Norme(IV_ActinidesComposition);
D = Distance( IV_ActinidesComposition / NormalisationFactor,
fFuelDataBank[i]);
if (D< distance)
{
distance = D;
N_BestEvolutionData = i;
}
}
DBGL
return fFuelDataBank[N_BestEvolutionData];
}
else if (fEvolutionDataInterpolation)
{
DBGL
map<double, int> distance_map = GetDistancesTo(isotopicvector, t);
map<double, int>::iterator it_distance;
int NClose = 64;
int Nstep = 0;
EvolutionData EvolInterpolate;
double SumOfDistance = 0;
for( it_distance = distance_map.begin(); it_distance != distance_map.end(); it_distance++)
{
double distance = (*it_distance).first;
int ED_Indice = (*it_distance).second;
if(distance == 0 )
{
EvolInterpolate = Multiply(1,fFuelDataBank[ED_Indice]);
return EvolInterpolate;
}
if(Nstep == 0)
EvolInterpolate = Multiply(1./distance, fFuelDataBank[ED_Indice]);
else
{
EvolutionData Evoltmp = EvolInterpolate;
EvolutionData Evoltmp2 = Multiply(1./distance, fFuelDataBank[ED_Indice]);
EvolInterpolate = Sum(Evoltmp, Evoltmp2);
Evoltmp.DeleteEvolutionData();
Evoltmp2.DeleteEvolutionData();
}
SumOfDistance += 1./distance;
Nstep++;
if(Nstep == NClose) break;
}
EvolutionData Evoltmp = EvolInterpolate;
EvolInterpolate.Clear();
EvolInterpolate = Multiply(1/SumOfDistance, Evoltmp);
Evoltmp.DeleteEvolutionData();
DBGL
return EvolInterpolate;
}
else
{
DBGL
IsotopicVector ActinidesCompositionAtT = fFuelDataBank[0].GetIsotopicVectorAt(t).GetActinidesComposition();
IsotopicVector IV_ActinidesComposition = isotopicvector.GetActinidesComposition();
double NormalisationFactor = Norme(IV_ActinidesComposition) / Norme( ActinidesCompositionAtT );
distance = Distance( IV_ActinidesComposition,
ActinidesCompositionAtT / NormalisationFactor,
fDistanceType,
fDistanceParameter);
for( int i = 1; i < (int)fFuelDataBank.size(); i++)
{
double D = 0;
ActinidesCompositionAtT = fFuelDataBank[i].GetIsotopicVectorAt(t).GetActinidesComposition();
IV_ActinidesComposition = isotopicvector.GetActinidesComposition();
NormalisationFactor = Norme(IV_ActinidesComposition) / Norme( ActinidesCompositionAtT );
D = Distance( IV_ActinidesComposition,
ActinidesCompositionAtT / NormalisationFactor,
fDistanceType,
fDistanceParameter);
if (D< distance)
{
distance = D;
N_BestEvolutionData = i;
}
}
INFO << distance << endl;
DBGL
return fFuelDataBank[N_BestEvolutionData];
}
}
//________________________________________________________________________
void XSM_CLOSEST::CalculateDistanceParameter()
{
DBGL
if(fDistanceType!=1){
WARNING << " Distance Parameter will be calculate even if the distance type is not the good one. Any Distance Parameters given by the user will be overwriten" << endl;
}
fDistanceParameter.Clear();
//We calculate the weight for the distance calculation.
int NevolutionDatainFuelDataBank = 0;
for( int i = 0; i < (int)fFuelDataBank.size(); i++)
{
NevolutionDatainFuelDataBank++;
map<ZAI ,double>::iterator itit;
map<ZAI ,double> isovector = fFuelDataBank[i].GetIsotopicVectorAt(0).GetIsotopicQuantity();
for(itit=isovector.begin(); itit != isovector.end(); itit++) //Boucle sur ZAI
{
double TmpXS=0;
for( int i=1; i<4; i++ ) //Loop on Reactions 1==fission, 2==capture, 3==n2n
TmpXS+= fFuelDataBank[i].GetXSForAt(0, (*itit).first, i);
fDistanceParameter.Add((*itit).first,TmpXS);
}
}
fDistanceParameter.Multiply( (double)1.0/NevolutionDatainFuelDataBank );
if(GetLog())
{
INFO <<"!!INFO!! Distance Parameters "<<endl;
map<ZAI ,double >::iterator it2;
for(it2 = fDistanceParameter.GetIsotopicQuantity().begin();it2 != fDistanceParameter.GetIsotopicQuantity().end(); it2++)
{
INFO << (*it2).first.Z() << " ";
INFO << (*it2).first.A() << " ";
INFO << (*it2).first.I() << " ";
INFO << ": " << (*it2).second;
INFO << endl;
}
INFO << endl;
}
DBGL
}
//________________________________________________________________________
void XSM_CLOSEST::SetDistanceParameter(IsotopicVector DistanceParameter)
{
fDistanceParameter = DistanceParameter;
INFO <<"!!INFO!! Distance Parameters "<<endl;
map<ZAI ,double >::iterator it2;
for(it2 = fDistanceParameter.GetIsotopicQuantity().begin();it2 != fDistanceParameter.GetIsotopicQuantity().end(); it2++)
{
INFO << (*it2).first.Z() << " ";
INFO << (*it2).first.A() << " ";
INFO << (*it2).first.I() << " ";
INFO << ": " << (*it2).second;
INFO << endl;
}
INFO << endl;
}
//________________________________________________________________________
void XSM_CLOSEST::SetDistanceType(int DistanceType)
{
fDistanceType=DistanceType;
if(fDistanceType==1){
CalculateDistanceParameter();
}
else if(fDistanceType == 2 && Norme(fDistanceParameter)==0){
// This is so bad!! You will probably unsynchronize all the reactor....
WARNING << " Distance use weight defined by user for each isotope, but no weight have been given" << endl<<"Use SetDistanceParameter()"<<endl;
exit(1);
}
else if (fDistanceType != 0 && fDistanceType != 1 && fDistanceType != 2 ){
ERROR << " Distancetype defined by the user isn't recognized by the code"<<endl;
exit(1);
}
}
\ No newline at end of file
source/branches/CLASSV3/Model/XS/XSM_CLOSEST.hxx deleted 100644 → 0
+ 0
− 134
View file @ 05226e12
#ifndef _XSM_CLOSEST_HXX
#define _XSM_CLOSEST_HXX
/*!
\file
\brief Header file for XSM_MLP_PWR_MOX class.
@authors BaM,BLG
@version 1.0
*/
#include "XSModel.hxx"
#include <string>
#include <fstream>
#include <iostream>
#include <map>
#include <vector>
typedef long long int cSecond;
using namespace std;
//-----------------------------------------------------------------------------//
/*!
Define a XSM_CLOSEST.
CLASS to get cross sections from a set of pre-calculation
(with MURE,or other depletion code)
get cross sections of the closest (in composition) calculation
@authors BaM,BLG
@version 1.0
*/
//________________________________________________________________________
class XSM_CLOSEST : public XSModel
{
public :
/*!
\name Constructor/Desctructor
*/
//@{
XSM_CLOSEST(string DB_index_file, bool oldreadmethod = false );
XSM_CLOSEST(CLASSLogger* Log, string DB_index_file, bool oldreadmethod = false );
~XSM_CLOSEST();
//{
//********* Get Method *********//
/*!
\name Get Method
*/
//@{
virtual EvolutionData GetCrossSections(IsotopicVector isotopicvector,double t=0) ; //!< Reason to live of this CLASS Return the closest Evolutiondata
vector< EvolutionData > GetFuelDataBank() const { return fFuelDataBank; } //!< Return the FuelDataBank
string GetDataBaseIndex() const { return fDataBaseIndex; } //!< Return the index Name
string GetFuelType() const { return fFuelType; } //!< Return the fuel type of the DB
vector<double> GetFuelParameter() const { return fFuelParameter; } //!< Return the Fuel parameter of the DB
pair<double,double> GetBurnUpRange() const { return fBurnUpRange;} //!< Return the BurnUp range of the DB
bool IsDefine(IsotopicVector IV) const; //!< True the key is define, false unstead
map<double, int> GetDistancesTo(IsotopicVector isotopicvector, double t = 0); //! Return a map containing the distance of each EvolutionData in the DataBase to the set IV at the t time
//@}
//********* Set Method *********//
/*!
\name Set Method
*/
//@{
void SetFuelDataBank(vector< EvolutionData > mymap) { fFuelDataBank = mymap; } //!< Set the FuelDataBank map
void SetDataBaseIndex(string database) { fDataBaseIndex = database;; ReadDataBase(); } //!< Set the Name of the database index
void SetOldReadMethod(bool val) { fOldReadMethod = val; ReadDataBase();} ///< use the old reading method
void SetWeightedDistanceCalculation(bool val = true) { fWeightedDistance = val;} ///< Set weighted Distance calculation
void SetInventoryEvolutionInterpolation(bool val = true) { fEvolutionDataInterpolation = val;} ///< Set weighted Distance calculation
void SetDistanceParameter(IsotopicVector DistanceParameter); ///< Define mannually the weight for each ZAI in the distance calculation
//{
/// Define the way to decide if two isotopic vectors are close.
/*!
// The different algorythm are:
// \li 0 is for the standard norme,
// \li 1 for each ZAI weighted with its XS,
// \li 2 for each ZAI weighted with coefficient given by the user.
*/
void SetDistanceType(int DistanceType);
//}
//********* Other Method *********//
/*!
\name Other Method
*/
//@{
void ReadDataBase(); ///< read the index file and fill the evolutionData map
void CalculateDistanceParameter(); ///< Calcul of the weight for each ZAI in the distance calculation from the mean XS of the FuelDataBank
//@}
private :
vector< EvolutionData > fFuelDataBank; ///< DataBanck map
string fDataBaseIndex; ///< Name of the index
bool fOldReadMethod; ///< use old DB format
bool fWeightedDistance; ///< USe XS weighted distance calculation
bool fEvolutionDataInterpolation; ///< USe XS weighted distance calculation
string fFuelType; ///< Type of fuel of the FuelDataBank
pair<double,double> fBurnUpRange; ///< Range of the Burn-up range of the FuelDataBank
vector<double> fFuelParameter; ///< Parameter needed by the equivalence model
int fDistanceType; ///< Set the distance calculation algorytm
/// \li 0 is for the standard norm (Default = 0),
/// \li 1 for each ZAI weighted with its XS,
/// \li 2 for each ZAI weighted with coefficient given by the user.
IsotopicVector fDistanceParameter; ///< weight for each ZAI in the distance calculation
};
#endif
source/branches/CLASSV3/Model/XS/XSM_MLP.cxx deleted 100644 → 0
+ 0
− 551
View file @ 05226e12
#include "XSModel.hxx"
#include "XSM_MLP.hxx"
#include "CLASSLogger.hxx"
#include "StringLine.hxx"
#include "TMVA/Reader.h"
#include "TMVA/Tools.h"
#include "TMVA/MethodCuts.h"
#include <TGraph.h>
#include <TString.h>
#include "TFile.h"
#include "TSystem.h"
#include "TROOT.h"
#include "TStopwatch.h"
#include <dirent.h>
#include <errno.h>
#include <sstream>
#include <string>
#include <iostream>
#include <fstream>
#include <algorithm>
#include <cmath>
//________________________________________________________________________
//
// XSM_MLP
//
//
//
//
//________________________________________________________________________
XSM_MLP::XSM_MLP(string TMVA_Weight_Directory,string InformationFile, bool IsTimeStep):XSModel(new CLASSLogger("XSM_MLP.log"))
{
fIsStepTime=IsTimeStep;
fTMVAWeightFolder = TMVA_Weight_Directory;
if(InformationFile=="")
fMLPInformationFile = TMVA_Weight_Directory+"/Data_Base_Info.nfo";
else
fMLPInformationFile=fTMVAWeightFolder+InformationFile;
GetMLPWeightFiles();
GetDataBaseInformation();
INFO<<"__A cross section interpolator using" <<endl;
INFO<<"Multi Layer Perceptron has been define__"<<endl;
INFO << " \t His TMVA folder is : \" " << fTMVAWeightFolder << "\"" << endl;
}
//________________________________________________________________________
XSM_MLP::XSM_MLP(CLASSLogger* Log,string TMVA_Weight_Directory,string InformationFile, bool IsTimeStep):XSModel(Log)
{
fIsStepTime=IsTimeStep;
fTMVAWeightFolder = TMVA_Weight_Directory;
if(InformationFile=="")
fMLPInformationFile = TMVA_Weight_Directory+"/Data_Base_Info.nfo";
else
fMLPInformationFile=fTMVAWeightFolder+InformationFile;
GetMLPWeightFiles();
GetDataBaseInformation();
INFO<<"__A cross section interpolator using" <<endl;
INFO<<"Multi Layer Perceptron has been define__"<<endl;
INFO << " \t His TMVA folder is : \" " << fTMVAWeightFolder << "\"" << endl;
}
//________________________________________________________________________
XSM_MLP::~XSM_MLP()
{
fMapOfTMVAVariableNames.clear();
}
//________________________________________________________________________
void XSM_MLP::GetDataBaseInformation()
{
ifstream FILE(fMLPInformationFile.c_str());
if(FILE.good())
{
while(!FILE.eof())
{
string line;
getline(FILE, line);
size_t foundRType = line.find("ReactorType :");
size_t foundFType = line.find("FuelType :");
size_t foundHM = line.find("Heavy Metal (t) :");
size_t foundPower = line.find("Thermal Power (W) :");
size_t foundTime = line.find("Time (s) :");
size_t foundZAI = line.find("Z A I Name (input MLP) :");
size_t foundDomain = line.find("Fuel range (Z A I min max) :");
int pos=0;
if(foundRType != std::string::npos)
{ StringLine::NextWord(line,pos,':');
fDataBaseRType = atof( (StringLine::NextWord(line,pos,':')).c_str() );
}
pos=0;
if(foundFType != std::string::npos)
{ StringLine::NextWord(line,pos,':');
fDataBaseFType = atof( (StringLine::NextWord(line,pos,':')).c_str() );
}
pos=0;
if(foundHM != std::string::npos)
{ StringLine::NextWord(line,pos,':');
fDataBaseHMMass = atof( (StringLine::NextWord(line,pos,':')).c_str() );
}
pos=0;
if(foundPower !=std::string::npos)
{ StringLine::NextWord(line,pos,':');
fDataBasePower = atof( (StringLine::NextWord(line,pos,':') ).c_str() );
}
pos=0;
if(foundTime!=std::string::npos)
{
StringLine::NextWord(line,pos,':');
while( pos< (int)line.size() )
fMLP_Time.push_back( atof( (StringLine::NextWord(line,pos,' ')).c_str() ));
}
pos=0;
if(foundZAI != std::string::npos)
{ string Z;
string A;
string I;
string Name;
int posoflinebeforbadline=0;
do
{ posoflinebeforbadline = FILE.tellg();
getline(FILE, line);
stringstream ssline;
ssline<<line;
ssline>>Z>>A>>I>>Name;
if(StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I) )
{
fMapOfTMVAVariableNames.insert( pair<ZAI,string>(ZAI(atoi(Z.c_str()),atoi(A.c_str()),atoi(I.c_str())),Name) );
}
}while((StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I)) && !FILE.eof());
FILE.seekg(posoflinebeforbadline); //return one line before
}
if(foundDomain != std::string::npos)
{ string Z;
string A;
string I;
string min;
string max;
int posoflinebeforbadline=0;
do
{ posoflinebeforbadline = FILE.tellg();
getline(FILE, line);
stringstream ssline;
ssline<<line;
ssline>>Z>>A>>I>>min>>max;
if(StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I) && StringLine::IsDouble(min) && StringLine::IsDouble(max) )
{
fFreshFuelDomain.insert( pair<ZAI,pair<double,double> >(ZAI(atoi(Z.c_str()),atoi(A.c_str()),atoi(I.c_str())),make_pair(atof(min.c_str()),atof(max.c_str()))) );
}
}
while((StringLine::IsDouble(Z) && StringLine::IsDouble(A) && StringLine::IsDouble(I) )&& !FILE.eof());
FILE.seekg(posoflinebeforbadline); //return one line before
}
}
}
else
{
ERROR << "Can't find/open file " << fMLPInformationFile << endl;
exit(0);
}
/********DEBUG*************************************/
INFO<<"\tMLP XS Data Base Information : "<<endl;
INFO<<"\t\tHeavy Metal (t) :"<<fDataBaseHMMass<<endl;
INFO<<"\t\tThermal Power (W) :"<<fDataBasePower<<endl;
INFO<<"\t\tTime (s) :"<<endl;
for (int i = 0; i < (int)fMLP_Time.size(); ++i)
INFO<<"\t\t\t"<<fMLP_Time[i]<<endl;
INFO<<"\t\tZ A I Name (input MLP) :"<<endl;
map<ZAI ,string >::iterator it;
for (it= fMapOfTMVAVariableNames.begin();it!=fMapOfTMVAVariableNames.end();it++)
INFO<<"\t\t\t"<< it->first.Z()<<" "<<it->first.A()<<" "<<it->second<<endl;
INFO<<"\t\tFuel range"<<endl;
for (map<ZAI,pair<double,double> >::iterator it_dom = fFreshFuelDomain.begin();it_dom!=fFreshFuelDomain.end();it_dom++)
INFO<<"\t\t\t"<< it_dom->second.first<<" <= "<<it_dom->first.Z()<<" "<<it_dom->first.A()<<" "<<it_dom->first.I()<<" <= "<<it_dom->second.second<<endl;;
}
//________________________________________________________________________
void XSM_MLP::GetMLPWeightFiles()
{DBGL
/**********Get All TMVA weight files*******************/
//check if the folder containing weights exists
DIR* rep = NULL;
struct dirent* fichierLu = NULL;
rep = opendir(fTMVAWeightFolder.c_str());
if (rep == NULL)
{
perror("");
exit(1);
}
/**Save file names of TMVA weights*/
fWeightFiles.resize(0);
while ((fichierLu = readdir(rep)) != NULL)
{
string FileName= fichierLu->d_name ;
if(FileName != "." && FileName != ".." )
{
if(FileName[FileName.size()-3]=='x' && FileName[FileName.size()-2]=='m' && FileName[FileName.size()-1]=='l' && FileName[0]!='.' )
fWeightFiles.push_back(FileName);
}
}
DBGL
}
//________________________________________________________________________
//________________________________________________________________________
//
// Time (MLP take time as parameter)
//________________________________________________________________________
//________________________________________________________________________
void XSM_MLP::ReadWeightFile(string Filename, int &Z, int &A, int &I, int &Reaction)
{
Z=-1;
A=-1;
I=-1;
Reaction=-1;
size_t found = Filename.find("XS");
string NameJOB;
NameJOB=Filename.substr(found);
int pos=0;
StringLine::NextWord(NameJOB, pos, '_');
Z = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
A = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
I = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
string sReaction = (StringLine::NextWord(NameJOB,pos,'_') ).c_str() ;
size_t foundext = sReaction.find(".weights.xml");
sReaction = sReaction.substr(0,foundext);
if(sReaction=="fis")
Reaction=0;
if(sReaction=="cap")
Reaction=1;
if(sReaction=="n2n")
Reaction=2;
if(Z<=0 || A<=0 || I<0 || Reaction==-1)
{
ERROR << " wrong TMVA weight format " << endl;
exit(0);
}
}
//________________________________________________________________________
TTree* XSM_MLP::CreateTMVAInputTree(IsotopicVector isotopicvector,int TimeStep)
{
/******Create Input data tree to be interpreted by TMVA::Reader***/
TTree* InputTree = new TTree("InTMP", "InTMP");
vector<float> InputTMVA;
for(int i = 0 ; i< (int)fMapOfTMVAVariableNames.size() ; i++)
InputTMVA.push_back(0);
float Time=0;
IsotopicVector IVInputTMVA;
map<ZAI ,string >::iterator it;
int j=0;
for( it = fMapOfTMVAVariableNames.begin() ; it!=fMapOfTMVAVariableNames.end() ; it++)
{
InputTree->Branch( ((*it).second).c_str() ,&InputTMVA[j], ((*it).second + "/F").c_str());
IVInputTMVA+= ((*it).first)*1;
j++;
}
if( !fIsStepTime)
InputTree->Branch( "Time" ,&Time ,"Time/F" );
IsotopicVector IVAccordingToUserInfoFile = isotopicvector.GetThisComposition(IVInputTMVA);
double Ntot = IVAccordingToUserInfoFile.GetSumOfAll();
IVAccordingToUserInfoFile = IVAccordingToUserInfoFile/Ntot;
j=0;
map<ZAI ,string >::iterator it2;
DBGV("INPUT TMVA");
for( it2 = fMapOfTMVAVariableNames.begin() ; it2!=fMapOfTMVAVariableNames.end() ; it2++)
{
InputTMVA[j] = IVAccordingToUserInfoFile.GetZAIIsotopicQuantity( (*it2).first ) ;
DBGV((*it2).first.Z()<<" "<<(*it2).first.A()<<" "<<InputTMVA[j]);
j++;
}
Time=fMLP_Time[TimeStep];
InputTree->Fill();
return InputTree;
}
//________________________________________________________________________
double XSM_MLP::ExecuteTMVA(string WeightFile,TTree* InputTree)
{
DBGV( "File :" << WeightFile);
// --- Create the Reader object
TMVA::Reader *reader = new TMVA::Reader( "Silent" );
// Create a set of variables and declare them to the reader
// - the variable names MUST corresponds in name and type to those given in the weight file(s) used
vector<float> InputTMVA;
for(int i = 0 ; i< (int)fMapOfTMVAVariableNames.size() ; i++)
InputTMVA.push_back(0);
Float_t Time;
map<ZAI ,string >::iterator it;
int j=0;
for( it = fMapOfTMVAVariableNames.begin() ; it!=fMapOfTMVAVariableNames.end() ; it++)
{ reader->AddVariable( ( (*it).second ).c_str(),&InputTMVA[j]);
j++;
}
if(!fIsStepTime)
reader->AddVariable( "Time" ,&Time);
// --- Book the MVA methods
string dir = fTMVAWeightFolder;
if(dir[dir.size()-1]!='/')
dir+="/";
// Book method MLP
TString methodName = "MLP method";
TString weightpath = dir + WeightFile ;
reader->BookMVA( methodName, weightpath );
map<ZAI ,string >::iterator it2;
j=0;
for( it2 = fMapOfTMVAVariableNames.begin() ; it2!=fMapOfTMVAVariableNames.end() ; it2++)
{
InputTree->SetBranchAddress(( (*it2).second ).c_str(),&InputTMVA[j]);
j++;
}
if(!fIsStepTime)
InputTree->SetBranchAddress( "Time" ,&Time );
InputTree->GetEntry(0);
Float_t val = (reader->EvaluateRegression( methodName ))[0];
delete reader;
DBGL
return (double)val;
}
//________________________________________________________________________
EvolutionData XSM_MLP::GetCrossSectionsTime(IsotopicVector IV)
{DBGL
EvolutionData EvolutionDataFromMLP = EvolutionData();
map<ZAI,TGraph*> ExtrapolatedXS[3];
/*************DATA BASE INFO****************/
EvolutionDataFromMLP.SetReactorType(fDataBaseRType);
EvolutionDataFromMLP.SetFuelType(fDataBaseFType);
EvolutionDataFromMLP.SetPower(fDataBasePower);
EvolutionDataFromMLP.SetHeavyMetalMass(fDataBaseHMMass);
/************* The Cross sections***********/
for(int i=0;i<int(fWeightFiles.size());i++)
{
int Z=-2;
int A=-2;
int I=-2;
int Reaction=-2;
ReadWeightFile( fWeightFiles[i], Z, A, I, Reaction);
if( Z >= GetZAIThreshold() )
{
for(int TimeStep=0;TimeStep<int(fMLP_Time.size());TimeStep++)
{
TTree* InputTree = CreateTMVAInputTree(IV,TimeStep);
pair< map<ZAI, TGraph*>::iterator, bool> IResult;
IResult = ExtrapolatedXS[Reaction].insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph() ) );
double XSValue = ExecuteTMVA(fWeightFiles[i],InputTree );
if(IResult.second )
{
(IResult.first)->second->SetPoint(0, (double)GetMLPTime()[TimeStep], XSValue );
}
else
{
(IResult.first)->second->SetPoint( (IResult.first)->second->GetN(), (double)GetMLPTime()[TimeStep], XSValue );
}
delete InputTree;
}
}
}
/**********Sorting TGraph*********/
for(int x=0;x<3;x++)
{ map<ZAI,TGraph*>::iterator it;
for(it = ExtrapolatedXS[x].begin(); it != ExtrapolatedXS[x].end(); it++)
it->second->Sort();
}
/**********Filling Matrices*/
EvolutionDataFromMLP.SetFissionXS(ExtrapolatedXS[0]);
EvolutionDataFromMLP.SetCaptureXS(ExtrapolatedXS[1]);
EvolutionDataFromMLP.Setn2nXS(ExtrapolatedXS[2]);
DBGL
return EvolutionDataFromMLP;
}
//________________________________________________________________________
//________________________________________________________________________
//
// Time step (1 MLP per time step)
//________________________________________________________________________
//________________________________________________________________________
void XSM_MLP::ReadWeightFileStep(string Filename, int &Z, int &A, int &I, int &Reaction, int &TimeStep)
{
Z=-1;
A=-1;
I=-1;
Reaction=-1;
size_t found = Filename.find("XS");
string NameJOB;
NameJOB=Filename.substr(found);
int pos=0;
StringLine::NextWord(NameJOB, pos, '_');
Z = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
A = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
I = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
string sReaction = (StringLine::NextWord(NameJOB,pos,'_') ).c_str() ;
if(sReaction=="fis")
Reaction=0;
if(sReaction=="cap")
Reaction=1;
if(sReaction=="n2n")
Reaction=2;
TimeStep = atof( (StringLine::NextWord(NameJOB,pos,'_') ).c_str() );
if(Z==-1 || A==-1 || I==-1 || Reaction==-1 || TimeStep==-1)
{
ERROR << " wrong TMVA weight format " << endl;
exit(0);
}
}
//________________________________________________________________________
EvolutionData XSM_MLP::GetCrossSectionsStep(IsotopicVector IV)
{DBGL
TTree* InputTree=CreateTMVAInputTree(IV);
//cout<<"=====Building Evolution Data From TMVA MLP====="<<endl;
EvolutionData EvolutionDataFromMLP = EvolutionData();
map<ZAI,TGraph*> ExtrapolatedXS[3];
/*************DATA BASE INFO****************/
EvolutionDataFromMLP.SetReactorType("PWR");
EvolutionDataFromMLP.SetFuelType("MOX");
EvolutionDataFromMLP.SetPower(fDataBasePower);
EvolutionDataFromMLP.SetHeavyMetalMass(fDataBaseHMMass);
/************* The Cross sections***********/
for(int i=0;i<int(fWeightFiles.size());i++)
{
int Z=-2;
int A=-2;
int I=-2;
int Reaction=-2;
int TimeStep=-2;
ReadWeightFileStep( fWeightFiles[i], Z, A, I, Reaction, TimeStep);
if( Z >= GetZAIThreshold() )
{
ZAI zaitmp = ZAI(Z,A,I);
pair< map<ZAI, TGraph*>::iterator, bool> IResult;
IResult = ExtrapolatedXS[Reaction].insert(pair<ZAI ,TGraph* >(ZAI(Z,A,I), new TGraph() ) );
if( IResult.second )
{
(IResult.first)->second->SetPoint(0, (double)GetMLPTime()[TimeStep], ExecuteTMVA(fWeightFiles[i],InputTree) );
}
else
{
(IResult.first)->second->SetPoint( (IResult.first)->second->GetN(), (double)GetMLPTime()[TimeStep], ExecuteTMVA(fWeightFiles[i],InputTree) );
}
}
}
/**********Sorting TGraph*********/
for(int x=0;x<3;x++)
{ map<ZAI,TGraph*>::iterator it;
for(it = ExtrapolatedXS[x].begin(); it != ExtrapolatedXS[x].end(); it++)
it->second->Sort();
}
/**********Filling Matrices*/
EvolutionDataFromMLP.SetFissionXS(ExtrapolatedXS[0]);
EvolutionDataFromMLP.SetCaptureXS(ExtrapolatedXS[1]);
EvolutionDataFromMLP.Setn2nXS(ExtrapolatedXS[2]);
delete InputTree;
DBGL
return EvolutionDataFromMLP;
}
//________________________________________________________________________
EvolutionData XSM_MLP::GetCrossSections(IsotopicVector IV ,double t)
{DBGL
if(t!=0)
WARNING << " Argument t has non effect here " << endl;
EvolutionData EV;
if(fIsStepTime)
EV=GetCrossSectionsStep(IV);
else
EV=GetCrossSectionsTime(IV);
DBGL
return EV;
}
//________________________________________________________________________
source/branches/CLASSV3/Model/XS/XSM_MLP.hxx deleted 100644 → 0
+ 0
− 91
View file @ 05226e12
#ifndef _XSM_MLP_HXX
#define _XSM_MLP_HXX
/*!
\file
\brief Header file for XSM_MLP class.
@authors BLG
@version 1.0
*/
#include "XSModel.hxx"
#include "TTree.h"
#include <string>
#include <fstream>
#include <iostream>
#include <map>
#include <vector>
typedef long long int cSecond;
using namespace std;
//-----------------------------------------------------------------------------//
/*!
Define a XSM_MLP.
This is the class to predict cross sections with a set of MultiLayerPerceptrons
@authors BLG
@version 1.0
*/
//________________________________________________________________________
class XSM_MLP : public XSModel
{
public :
/*!
\name Constructor/Desctructor
*/
//@{
/// CONSTRUCTOR
/*!
\param string TMVA_Weight_Directory The directory where all the TMVA weight are located
\param string InformationFile, Name of the information file located in TMVA_Weight_Directory (default: Data_Base_Info.nfo)
\param bool 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=true);
XSM_MLP(CLASSLogger* Log,string TMVA_Weight_Directory,string InformationFile="",bool IsTimeStep=false);
~XSM_MLP();
//{
virtual EvolutionData GetCrossSections(IsotopicVector IV,double t=0) ;//!< Return calculated cross section by the MLP regression
private :
void GetDataBaseInformation(); //<! Read information file and fill Reactor Type, Fuel type, HM mass, Power, time vector, and TMVA input variables names (looks at Data Bases example for format details)
vector<double> GetMLPTime() {return fMLP_Time; }//<! Return time vector (seconds) defined in the DataBaseInformation file
void GetMLPWeightFiles();//<! Find all .xml file in TMVA_Weight_Directory
void ReadWeightFile(string Filename, int &Z, int &A, int &I, int &Reaction) ;//<! Select the reaction according to the weight file name (file name is formated !!) read the manual to know it (formated for fIsTimeStep==false)
double ExecuteTMVA(string WeightFile, TTree* InputTree);//!<Execute the MLP according to the input tree created by CreateTMVAInputTree
TTree* CreateTMVAInputTree(IsotopicVector isotopicvector,int TimeStep=0);//!<Create input tmva tree to be read by ExecuteTMVA
EvolutionData GetCrossSectionsTime(IsotopicVector IV) ;//!< Return calculated cross section by the MLP regression when fIsTimeStep==false
void ReadWeightFileStep(string Filename, int &Z, int &A, int &I, int &Reaction, int &TimeStep) ;//<! Select the reaction according to the weight file name (file name is formated !!) read the manual to know it (formated for fIsTimeStep==true)
EvolutionData GetCrossSectionsStep(IsotopicVector IV) ;//!< Return calculated cross section by the MLP regression when fIsTimeStep==true
vector<double> fMLP_Time;//<! Time vector of the data base
vector<string> fWeightFiles;//<! All the weight file contains in fTMVAWeightFolder
string fTMVAWeightFolder;//<! folder containing all the weight file
string fMLPInformationFile;//<! file containing Reactor Type, Fuel type, HM mass, Power, time vector, and TMVA input variables names (looks the manual for format details)
double fDataBasePower;//<!Power of the data base (read from fMLPInformationFile )
double fDataBaseHMMass;//<!Heavy metal mass of the data base (read from fMLPInformationFile )
string fDataBaseFType; //<! Reactor Type (e.g PWR, RNR, ADS..)
string fDataBaseRType; //<! Fuel Type (e.g MOX, UOX, ThU, ThPu ...)
bool fIsStepTime;//<!true if one TMVA weihgt per step time is requiered otherwise it assumes time is part of the MLP inputs
map<ZAI,string> fMapOfTMVAVariableNames;//<! List of TMVA input variable names (read from fMLPInformationFile ) , name depends on the training step
};
#endif
source/branches/CLASSV3/include/CLASSBackEnd.hxx deleted 100644 → 0
+ 0
− 122
View file @ 05226e12
#ifndef _CLASSBACKEND_HXX
#define _CLASSBACKEND_HXX
/*!
\file
\brief Header file for CLASSFacility class.
*/
#include <string>
#include <fstream>
#include <map>
#include "CLASSFacility.hxx"
#include "IsotopicVector.hxx"
#include "DecayDataBank.hxx"
#include "TNamed.h"
using namespace std;
typedef long long int cSecond;
//-----------------------------------------------------------------------------//
/*!
Define a CLASS Facility.
The aim of these class is synthetyse all the commum properties of the nuclear facilities which are involve in the BackEnd Fuel cycle.
@author BaM
@version 2.0
*/
//________________________________________________________________________
class CLASSBackEnd : public CLASSFacility
{
public :
///< Normal Constructor.
CLASSBackEnd(int type = 0);
CLASSBackEnd(CLASSLogger* log,int type = 0);
CLASSBackEnd(CLASSLogger* log, cSecond cycletime, int type = 0);
//********* Get Method *********//
/*!
\name Get Function
*/
//@{
vector<IsotopicVector> GetIVArray() const { return fIVArray; } //!< Return the IsotopicVector Array
vector<cSecond> GetIVArrayArrivalTime() const { return fIVArrayArrivalTime;} //!<Return the pointer to the OUtBackEndFacility
int GetIVNumber() const { return fIVArray.size();}
bool GetStorageType() const { return fIsStorageType;} //!< Return the storageType
IsotopicVector GetIV(int i) const { if(i < (int)fIVArray.size()) return fIVArray[i];
else return IsotopicVector(); }
#ifndef __CINT__
DecayDataBank* GetDecayDataBank() { return fDecayDataBase;} //!< Return the pointer to the decay DataBank
CLASSBackEnd* GetOutBackEndFacility() const { return fOutBackEndFacility;} //!<Return the pointer to the OUtBackEndFacility
virtual map<cSecond,int> GetTheBackEndTimePath();
#endif
//@}
//********* Set Method *********//
/*!
\name Set Function
*/
//@{
void SetIsStorageType(bool val = true) { fIsStorageType = val;} //! Set the fIsStorage bool
virtual void SetIVArray(vector<IsotopicVector> ivarray) { fIVArray = ivarray; } //!< Set The isotopicVector Array
#ifndef __CINT__
void SetDecayDataBank(DecayDataBank* decayDB) { fDecayDataBase = decayDB;} //! Set the Decay DataBank
virtual void SetOutBackEndFacility(CLASSBackEnd* befacility) { fOutBackEndFacility = befacility;
fIsStorageType = false; } //! Set a Out Facility for the fuel
#endif
using CLASSFacility::SetName;
//@}
/*!
\name BackEndFacility specific Method
*/
//@{
virtual void AddIV(IsotopicVector isotopicvector); //!< Add an Isotopicvector to the IVArray
void ClearIVArray(); //!< Empty the IVArray removing all fuel stored
//@}
virtual void Evolution(cSecond t) {} //!< Performe the Evolution to the Time t
void UpdateInsideIV();
protected :
IsotopicVector GetDecay(IsotopicVector isotopicvector, cSecond t); //!< Get IsotopicVector Decay at the t time
vector<IsotopicVector> fIVArray; ///< Vector containning all the fuel stored.
vector<cSecond> fIVArrayArrivalTime; ///< Vector containning all the fuel stored.
#ifndef __CINT__
CLASSBackEnd* fOutBackEndFacility; //!< Facility getting the fuel at the end of the cycle
#endif
//********* Internal Parameter *********//
private :
bool fIsStorageType; //!< True if there is not OutBAckEndFacility (like a storage...)
#ifndef __CINT__
DecayDataBank* fDecayDataBase; //!< Pointer to the Decay DataBase
#endif
ClassDef(CLASSBackEnd,2);
};
#endif
source/branches/CLASSV3/include/CLASSConstante.hxx deleted 100644 → 0
+ 0
− 11
View file @ 05226e12
#ifndef _CLASSConstante_HXX_
#define _CLASSConstante_HXX_
//CLASS library
#include "ZAIMass.hxx"
const double AVOGADRO = 6.02214129e23;
const ZAIMass cZAIMass;
#endif
source/branches/CLASSV3/include/CLASSFacility.hxx deleted 100644 → 0
+ 0
− 141
View file @ 05226e12
#ifndef _CLASSFACILITY_HXX
#define _CLASSFACILITY_HXX
/*!
\file
\brief Header file for CLASSFacility class.
*/
#include <string>
#include <fstream>
#include "CLASSObject.hxx"
#include "IsotopicVector.hxx"
#include "DecayDataBank.hxx"
#include "TNamed.h"
using namespace std;
typedef long long int cSecond;
class Scenario;
//-----------------------------------------------------------------------------//
/*!
Define a CLASS Facility.
The aim of these class is to regroup all the commum properties of the nuclear facilities.
@author BaM
@version 2.0
*/
//________________________________________________________________________
class CLASSFacility : public CLASSObject
{
public :
///< Normal Constructor.
CLASSFacility(int type = 0);
CLASSFacility(CLASSLogger* log, int type = 0);
CLASSFacility(CLASSLogger* log, cSecond cycletime, int type = 0);
CLASSFacility(CLASSLogger* log, cSecond creationtime, cSecond lifetime, int type = 0);
CLASSFacility(CLASSLogger* log, cSecond startingtime, cSecond lifetime, cSecond cycletime, int type = 0);
//********* Get Method *********//
/*!
\name Get Function
*/
//@{
int GetId() const { return fId; } //!< Return the Facility Parc'Is
IsotopicVector GetInsideIV() const { return fInsideIV; } //!< Return the IV contain in the Facility
int GetFacilityType() const { return fFacilityType; }
cSecond GetInternalTime() const { return fInternalTime; } //!< Return Creation Time
cSecond GetCycleTime() const { return fCycleTime; } //!< Return the cycle time of the Facility
cSecond GetCreationTime() const { return fCreationTime; } //!< Return the creation time of the Facility
cSecond GetLifeTime() const { return fLifeTime; } //!< Return the life time of the Facility
Scenario* GetParc() { return fParc; } //!< return the pointer to the Park
IsotopicVector GetCumulativeIVIn() { return fCumulativeIVIn;} //!< return the culative sum of all incoming IV
IsotopicVector GetCumulativeIVOut() { return fCumulativeIVOut;} //!< return the culative sum of all outcoming IV
//@}
//********* Set Method *********//
/*!
\name Set Function
*/
//@{
void SetId(int id) { fId = id; } //!< Set The Facility Parc'Id
void SetParc(Scenario* parc) { fParc = parc; } //!< Set the Pointer to the Parc
void SetFacilityType(int type) { fFacilityType = type; } //!< Set the facility type :
/// \li 2 reactor Studown
/// \li 4 start/End of reactor cycle,
/// \li 8 end of Cooling,
/// \li 16 fuel Fabrication
using CLASSObject::SetName;
using CLASSObject::GetName;
void SetInsideIV(IsotopicVector isotopicvector) { fInsideIV = isotopicvector; } //!< Set the IV inside the Facility Core
void SetCreationTime(double creationtime) { fCreationTime = (cSecond)creationtime;} //!< Set the creation Time
void SetLifeTime(double lifetime) { fLifeTime = (cSecond)lifetime; } //!< Set the life time of the facility
virtual void SetCycleTime(double cycletime) { fCycleTime = (cSecond)cycletime; } //!< Set the cycle time (Cycle of the loading Plan)
void SetInCycleTime(double incycletime) { fInCycleTime = (cSecond)incycletime; fIsStarted = true; } //!< Set the cycle time (Cycle of the loading Plan)
void SetInternalTime(double internaltime) { fInternalTime = (cSecond)internaltime; } //!< Set the cycle time (Cycle of the loading Plan)
//@}
/*!
\name Evolution Method
*/
//@{
void AddCumulativeIVIn(IsotopicVector IV) { fCumulativeIVIn += IV;} //!< Add the Input IsotopicVector the The cumulative IV IN
void AddCumulativeIVOut(IsotopicVector IV) { fCumulativeIVOut += IV;} //!< Add the Input IsotopicVector the The cumulative IV OUT
virtual void Evolution(cSecond t) = 0; //!< Performe the Evolution to the Time t
virtual void Dump() { } //!< Write Modification (IV In/Out, filling the TF...)
//@}
protected :
bool fIsStarted; ///< True if Running, False Otherwise
bool fIsShutDown; ///< True if the facility is stoped, False Otherwise
bool fIsAtEndOfCycle; ///< True if Reaching the End of a Facility Cycle
cSecond fInternalTime; ///< Internal Clock
cSecond fInCycleTime; ///< Time spend since the beginning of the last Cycle
cSecond fCycleTime; ///< Cycle duration Time
IsotopicVector fInsideIV; ///< All IV in the Facility (fuel for reactor, total for all others...)
IsotopicVector fCumulativeIVIn; ///< All IV in the Facility (fuel for reactor, total for all others...)
IsotopicVector fCumulativeIVOut; ///< All IV in the Facility (fuel for reactor, total for all others...)
//********* Internal Parameter *********//
private :
int fId; //!< Identity of the Facility inside the Parc
int fFacilityType; ///< Type of facility :
/// \li 4 reactor,
/// \li 8 Pool,
/// \li 16 FabricationPlant.
Scenario* fParc; //!< Pointer to the main Parc
cSecond fCreationTime; ///< CLASS Universal Time of Creation
cSecond fLifeTime; ///< Time of life Of the Reactor (Operating's Duration)
ClassDef(CLASSFacility,1);
};
#endif
source/branches/CLASSV3/include/CLASSFuel.hxx deleted 100644 → 0
+ 0
− 62
View file @ 05226e12
#ifndef _CLASSFUEL_HXX
#define _CLASSFUEL_HXX
/*!
\file
\brief Header file for CLASSFuel class.
@author BaM
@version 2.0
*/
#include <string>
#include <fstream>
#include "CLASSObject.hxx"
#include "PhysicsModels.hxx"
#include "EvolutionData.hxx"
using namespace std;
//-----------------------------------------------------------------------------//
/*!
Define a CLASS Object.
The aim of these class is synthetyse all the commum properties to all CLASS Fuel Element.
@author BaM
@version 2.0
*/
//________________________________________________________________________
class CLASSFuel : public CLASSObject
{
public :
///< Normal Constructor.
CLASSFuel(EvolutionData* evo);
CLASSFuel(PhysicsModels* evo);
virtual CLASSFuel* Clone() { return new CLASSFuel(*this); } //!< Correct way to copy a CLASSFuel in case of derivation
EvolutionData* GetEvolutionData() {return fEvolutionData;}
PhysicsModels* GetPhysicsModels() {return fPhysicsModels;}
using CLASSObject::SetName;
using CLASSObject::GetName;
protected :
private :
EvolutionData* fEvolutionData;
PhysicsModels* fPhysicsModels;
};
#endif
source/branches/CLASSV3/include/CLASSFuelPlan.hxx deleted 100644 → 0
+ 0
− 65
View file @ 05226e12
#ifndef _CLASSFUELPLAN_HXX
#define _CLASSFUELPLAN_HXX
/*!
\file
\brief Header file for CLASSFuelPlan class.
@author BaM
@version 2.0
*/
#include <string>
#include <fstream>
#include <map>
#include "CLASSObject.hxx"
#include "CLASSFuel.hxx"
using namespace std;
typedef long long int cSecond;
//-----------------------------------------------------------------------------//
/*!
Define a CLASS Object.
The aim of these class is synthetyse all the commum properties to all CLASS Fuel Element.
@author BaM
@version 2.0
*/
//________________________________________________________________________
class CLASSFuelPlan : public CLASSObject
{
public :
///< Normal Constructor.
CLASSFuelPlan();
CLASSFuelPlan(CLASSLogger* log);
void AddFuel(cSecond time, CLASSFuel fuel, double BurnUp);
void AddFuel(cSecond time, EvolutionData* fuel, double BurnUp) {AddFuel( time, CLASSFuel(fuel), BurnUp);}
void AddFuel(cSecond time, PhysicsModels* fuel, double BurnUp) {AddFuel( time, CLASSFuel(fuel), BurnUp);}
pair< CLASSFuel, double> GetFuelAt(cSecond t);
using CLASSObject::SetName;
using CLASSObject::GetName;
protected :
private :
map< cSecond, pair< CLASSFuel, double > > fLoadingPlan; ///< Loading PLan to change the EvolutionData (and the associetedBurnup) according to the Plan
};
#endif
source/branches/CLASSV3/include/CLASSHeaders.hxx deleted 100755 → 0
+ 0
− 25
View file @ 05226e12
#ifndef _CLASSHEADERS_HXX_
#define _CLASSHEADERS_HXX_
//CLASS library
#include "CLASSConstante.hxx"
#include "Scenario.hxx"
#include "Reactor.hxx"
#include "Pool.hxx"
#include "FabricationPlant.hxx"
#include "SeparationPlant.hxx"
#include "Storage.hxx"
#include "IsotopicVector.hxx"
#include "ZAI.hxx"
#include "CLASSLogger.hxx"
#include "EvolutionData.hxx"
#include "PhysicsModels.hxx"
#include "DecayDataBank.hxx"
#include "ZAIMass.hxx"
#endif
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment