diff --git a/source/trunk/Model/Irradiation/IM_Matrix_Decay.cxx b/source/trunk/Model/Irradiation/IM_Matrix_Decay.cxx index a74c953ae1bfe5de9d60f643af21830c6a74e603..fe7058ce619e2a93e374b9505080b8bb1159aef5 100644 --- a/source/trunk/Model/Irradiation/IM_Matrix_Decay.cxx +++ b/source/trunk/Model/Irradiation/IM_Matrix_Decay.cxx @@ -15,7 +15,7 @@ #include <TGraph.h> #include <TString.h> - +#include <TMatrixT.h> #include <sstream> #include <string> @@ -93,13 +93,21 @@ TMatrixT<double> IM_Matrix_Decay::ExponentialCalculation(TMatrixT<double> myMatr } -IsotopicVector IM_Matrix_Decay::GetDecay(IsotopicVector Mother_IV, time) +IsotopicVector IM_Matrix_Decay::GetDecay(IsotopicVector Mother_IV, double time) { + IsotopicVector DecayIV; - TMatrixT<double> NMatrix = TMatrixT<double>(decayindex.size(),1)) + TMatrixT<double> NMatrix_0 = TMatrixT<double>(fReverseMatrixIndex.size(),1); + TMatrixT<double> NMatrix_t = TMatrixT<double>(fReverseMatrixIndex.size(),1); for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++) - NMatrix[i] = 0; + { + NMatrix_0[i] = 0; + NMatrix_t[i] = 0; + } + + + { // Filling the t=0 State; map<ZAI, double > isotopicquantity = Mother_IV.GetIsotopicQuantity(); @@ -119,7 +127,7 @@ IsotopicVector IM_Matrix_Decay::GetDecay(IsotopicVector Mother_IV, time) it2 = fMatrixIndex.find( ZAI(-3,-3,-3) ); - NMatrix[ (*it2).second ][0] = (*it).second ; + NMatrix_0[ (*it2).second ][0] = (*it).second ; } @@ -128,6 +136,23 @@ IsotopicVector IM_Matrix_Decay::GetDecay(IsotopicVector Mother_IV, time) } + TMatrixT<double> exp_T_Matrix = 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) exp_T_Matrix[j][k] = exp(time); + else exp_T_Matrix[j][k] = 0; + } + + NMatrix_t = fExponentialDecayMatrix * exp_T_Matrix * NMatrix_0; + + for(int i = 0; i < (int)fReverseMatrixIndex.size(); i++) + { + DecayIV += fReverseMatrixIndex[i]* NMatrix_t[i][0]; + + } + return DecayIV; + } diff --git a/source/trunk/Model/Irradiation/IM_Matrix_Decay.hxx b/source/trunk/Model/Irradiation/IM_Matrix_Decay.hxx index f38ee66948324e1b78ea957ec308d36a60ac3b1a..2931d856b1f1d12d8f208d1d238a76e919c1e396 100644 --- a/source/trunk/Model/Irradiation/IM_Matrix_Decay.hxx +++ b/source/trunk/Model/Irradiation/IM_Matrix_Decay.hxx @@ -62,7 +62,7 @@ class IM_Matrix_Decay : public IrradiationModel //@} - IsotopicVector GetDecay(IsotopicVector Mother_IV, time); + IsotopicVector GetDecay(IsotopicVector Mother_IV, double time); TMatrixT<double> ExponentialCalculation(TMatrixT<double> myMatrix);