From 540ed361719893594c89ca6d9c95f4db78b80aa1 Mon Sep 17 00:00:00 2001 From: "Bastian M. Wojek" Date: Tue, 16 Mar 2010 20:44:57 +0000 Subject: [PATCH] Added lambda(T) and 1/lambda(T) operators to libGapIntegrals in order to increase flexibility --- .../libGapIntegrals/TGapIntegrals.cpp | 320 +++++++++++++++++- src/external/libGapIntegrals/TGapIntegrals.h | 184 ++++++++++ .../libGapIntegrals/TGapIntegralsLinkDef.h | 13 + 3 files changed, 510 insertions(+), 7 deletions(-) diff --git a/src/external/libGapIntegrals/TGapIntegrals.cpp b/src/external/libGapIntegrals/TGapIntegrals.cpp index f27b6482..092047a2 100644 --- a/src/external/libGapIntegrals/TGapIntegrals.cpp +++ b/src/external/libGapIntegrals/TGapIntegrals.cpp @@ -48,13 +48,28 @@ ClassImp(TGapNonMonDWave1) ClassImp(TGapNonMonDWave2) ClassImp(TGapPowerLaw) +ClassImp(TLambdaSWave) +ClassImp(TLambdaDWave) +ClassImp(TLambdaAnSWave) +ClassImp(TLambdaNonMonDWave1) +ClassImp(TLambdaNonMonDWave2) +ClassImp(TLambdaPowerLaw) + +ClassImp(TLambdaInvSWave) +ClassImp(TLambdaInvDWave) +ClassImp(TLambdaInvAnSWave) +ClassImp(TLambdaInvNonMonDWave1) +ClassImp(TLambdaInvNonMonDWave2) +ClassImp(TLambdaInvPowerLaw) + +ClassImp(TFilmMagnetizationDWave) + // s wave gap integral TGapSWave::TGapSWave() { TGapIntegral *gapint = new TGapIntegral(); fGapIntegral = gapint; gapint = 0; - delete gapint; fTemp.clear(); fTempIter = fTemp.end(); @@ -67,7 +82,6 @@ TGapDWave::TGapDWave() { TDWaveGapIntegralCuhre *gapint = new TDWaveGapIntegralCuhre(); fGapIntegral = gapint; gapint = 0; - delete gapint; fTemp.clear(); fTempIter = fTemp.end(); @@ -80,7 +94,6 @@ TGapAnSWave::TGapAnSWave() { TAnSWaveGapIntegralCuhre *gapint = new TAnSWaveGapIntegralCuhre(); fGapIntegral = gapint; gapint = 0; - delete gapint; fTemp.clear(); fTempIter = fTemp.end(); @@ -93,7 +106,6 @@ TGapNonMonDWave1::TGapNonMonDWave1() { TNonMonDWave1GapIntegralCuhre *gapint = new TNonMonDWave1GapIntegralCuhre(); fGapIntegral = gapint; gapint = 0; - delete gapint; fTemp.clear(); fTempIter = fTemp.end(); @@ -106,7 +118,6 @@ TGapNonMonDWave2::TGapNonMonDWave2() { TNonMonDWave2GapIntegralCuhre *gapint = new TNonMonDWave2GapIntegralCuhre(); fGapIntegral = gapint; gapint = 0; - delete gapint; fTemp.clear(); fTempIter = fTemp.end(); @@ -115,10 +126,50 @@ TGapNonMonDWave2::TGapNonMonDWave2() { fPar.clear(); } +TLambdaSWave::TLambdaSWave() { + fLambdaInvSq = new TGapSWave(); +} + +TLambdaDWave::TLambdaDWave() { + fLambdaInvSq = new TGapDWave(); +} + +TLambdaAnSWave::TLambdaAnSWave() { + fLambdaInvSq = new TGapAnSWave(); +} + +TLambdaNonMonDWave1::TLambdaNonMonDWave1() { + fLambdaInvSq = new TGapNonMonDWave1(); +} + +TLambdaNonMonDWave2::TLambdaNonMonDWave2() { + fLambdaInvSq = new TGapNonMonDWave2(); +} + +TLambdaInvSWave::TLambdaInvSWave() { + fLambdaInvSq = new TGapSWave(); +} + +TLambdaInvDWave::TLambdaInvDWave() { + fLambdaInvSq = new TGapDWave(); +} + +TLambdaInvAnSWave::TLambdaInvAnSWave() { + fLambdaInvSq = new TGapAnSWave(); +} + +TLambdaInvNonMonDWave1::TLambdaInvNonMonDWave1() { + fLambdaInvSq = new TGapNonMonDWave1(); +} + +TLambdaInvNonMonDWave2::TLambdaInvNonMonDWave2() { + fLambdaInvSq = new TGapNonMonDWave2(); +} + TGapSWave::~TGapSWave() { delete fGapIntegral; fGapIntegral = 0; - + fTemp.clear(); fTempIter = fTemp.end(); fIntegralValues.clear(); @@ -170,6 +221,57 @@ TGapNonMonDWave2::~TGapNonMonDWave2() { fPar.clear(); } +TLambdaSWave::~TLambdaSWave() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaDWave::~TLambdaDWave() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaAnSWave::~TLambdaAnSWave() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaNonMonDWave1::~TLambdaNonMonDWave1() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaNonMonDWave2::~TLambdaNonMonDWave2() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaInvSWave::~TLambdaInvSWave() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaInvDWave::~TLambdaInvDWave() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaInvAnSWave::~TLambdaInvAnSWave() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaInvNonMonDWave1::~TLambdaInvNonMonDWave1() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + +TLambdaInvNonMonDWave2::~TLambdaInvNonMonDWave2() { + delete fLambdaInvSq; + fLambdaInvSq = 0; +} + + double TGapSWave::operator()(double t, const vector &par) const { assert(par.size() == 2); // two parameters: Tc, Delta(0) @@ -309,7 +411,7 @@ double TGapDWave::operator()(double t, const vector &par) const { double TGapAnSWave::operator()(double t, const vector &par) const { - assert(par.size() == 3); // two parameters: Tc, Delta(0), a + assert(par.size() == 3); // three parameters: Tc, Delta(0), a if (t<=0.0) return 1.0; @@ -529,3 +631,207 @@ double TGapPowerLaw::operator()(double t, const vector &par) const { return 1.0 - pow(t/par[0], par[1]); } + + +double TLambdaSWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 2); // two parameters: Tc, Delta0 + + if (t >= par[0]) + return -1.0; + + if (t <= 0.0) + return 1.0; + + return 1.0/sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaDWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 2); // two parameters: Tc, Delta0 + + if (t >= par[0]) + return -1.0; + + if (t <= 0.0) + return 1.0; + + return 1.0/sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaAnSWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 3); // three parameters: Tc, Delta0, a + + if (t >= par[0]) + return -1.0; + + if (t <= 0.0) + return 1.0; + + return 1.0/sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaNonMonDWave1::operator()(double t, const vector &par) const +{ + assert(par.size() == 3); // three parameters: Tc, Delta0, a + + if (t >= par[0]) + return -1.0; + + if (t <= 0.0) + return 1.0; + + return 1.0/sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaNonMonDWave2::operator()(double t, const vector &par) const +{ + assert(par.size() == 3); // three parameters: Tc, Delta0, a + + if (t >= par[0]) + return -1.0; + + if (t <= 0.0) + return 1.0; + + return 1.0/sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaPowerLaw::operator()(double t, const vector &par) const { + + assert(par.size() == 2); // two parameters: Tc, N + + if(t <= 0.0) + return 1.0; + else if (t >= par[0]) + return -1.0; + + return 1.0/sqrt(1.0 - pow(t/par[0], par[1])); + +} + + +double TLambdaInvSWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 2); // two parameters: Tc, Delta0 + + if (t >= par[0]) + return 0.0; + + if (t <= 0.0) + return 1.0; + + return sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaInvDWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 2); // two parameters: Tc, Delta0 + + if (t >= par[0]) + return 0.0; + + if (t <= 0.0) + return 1.0; + + return sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaInvAnSWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 3); // three parameters: Tc, Delta0, a + + if (t >= par[0]) + return 0.0; + + if (t <= 0.0) + return 1.0; + + return sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaInvNonMonDWave1::operator()(double t, const vector &par) const +{ + assert(par.size() == 3); // three parameters: Tc, Delta0, a + + if (t >= par[0]) + return 0.0; + + if (t <= 0.0) + return 1.0; + + return sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaInvNonMonDWave2::operator()(double t, const vector &par) const +{ + assert(par.size() == 3); // three parameters: Tc, Delta0, a + + if (t >= par[0]) + return 0.0; + + if (t <= 0.0) + return 1.0; + + return sqrt((*fLambdaInvSq)(t, par)); + +} + +double TLambdaInvPowerLaw::operator()(double t, const vector &par) const { + + assert(par.size() == 2); // two parameters: Tc, N + + if(t <= 0.0) + return 1.0; + else if (t >= par[0]) + return 0.0; + + return sqrt(1.0 - pow(t/par[0], par[1])); + +} + +TFilmMagnetizationDWave::TFilmMagnetizationDWave() +{ + fLambdaInvSq = new TGapDWave(); + fPar.clear(); +} + +TFilmMagnetizationDWave::~TFilmMagnetizationDWave() +{ + delete fLambdaInvSq; + fLambdaInvSq = 0; + fPar.clear(); +} + +double TFilmMagnetizationDWave::operator()(double t, const vector &par) const +{ + assert(par.size() == 4); // four parameters: Tc, Delta0, lambda0, film-thickness + + fPar = par; + + if (t >= fPar[0]) + return 0.0; + + vector parForGapIntegral; + parForGapIntegral.push_back(fPar[0]); + parForGapIntegral.push_back(fPar[1]); + + double d_2l(0.5*fPar[3]/fPar[2]*sqrt((*fLambdaInvSq)(t, parForGapIntegral))); + + parForGapIntegral.clear(); + + return tanh(d_2l)/d_2l - 1.0; + +} + + diff --git a/src/external/libGapIntegrals/TGapIntegrals.h b/src/external/libGapIntegrals/TGapIntegrals.h index a415b826..3c4c54e9 100644 --- a/src/external/libGapIntegrals/TGapIntegrals.h +++ b/src/external/libGapIntegrals/TGapIntegrals.h @@ -154,5 +154,189 @@ private: ClassDef(TGapPowerLaw,1) }; +class TLambdaSWave : public PUserFcnBase { + +public: + TLambdaSWave(); + ~TLambdaSWave(); + + double operator()(double, const vector&) const; + +private: + TGapSWave *fLambdaInvSq; + + ClassDef(TLambdaSWave,1) +}; + +class TLambdaDWave : public PUserFcnBase { + +public: + TLambdaDWave(); + ~TLambdaDWave(); + + double operator()(double, const vector&) const; + +private: + TGapDWave *fLambdaInvSq; + + ClassDef(TLambdaDWave,1) +}; + +class TLambdaAnSWave : public PUserFcnBase { + +public: + TLambdaAnSWave(); + ~TLambdaAnSWave(); + + double operator()(double, const vector&) const; + +private: + TGapAnSWave *fLambdaInvSq; + + ClassDef(TLambdaAnSWave,1) +}; + +class TLambdaNonMonDWave1 : public PUserFcnBase { + +public: + TLambdaNonMonDWave1(); + ~TLambdaNonMonDWave1(); + + double operator()(double, const vector&) const; + +private: + TGapNonMonDWave1 *fLambdaInvSq; + + ClassDef(TLambdaNonMonDWave1,1) +}; + +class TLambdaNonMonDWave2 : public PUserFcnBase { + +public: + TLambdaNonMonDWave2(); + ~TLambdaNonMonDWave2(); + + double operator()(double, const vector&) const; + +private: + TGapNonMonDWave2 *fLambdaInvSq; + + ClassDef(TLambdaNonMonDWave2,1) +}; + + +class TLambdaPowerLaw : public PUserFcnBase { + +public: + TLambdaPowerLaw() {} + ~TLambdaPowerLaw() {} + + double operator()(double, const vector&) const; + +private: + + ClassDef(TLambdaPowerLaw,1) +}; + +class TLambdaInvSWave : public PUserFcnBase { + +public: + TLambdaInvSWave(); + ~TLambdaInvSWave(); + + double operator()(double, const vector&) const; + +private: + TGapSWave *fLambdaInvSq; + + ClassDef(TLambdaInvSWave,1) +}; + +class TLambdaInvDWave : public PUserFcnBase { + +public: + TLambdaInvDWave(); + ~TLambdaInvDWave(); + + double operator()(double, const vector&) const; + +private: + TGapDWave *fLambdaInvSq; + + ClassDef(TLambdaInvDWave,1) +}; + +class TLambdaInvAnSWave : public PUserFcnBase { + +public: + TLambdaInvAnSWave(); + ~TLambdaInvAnSWave(); + + double operator()(double, const vector&) const; + +private: + TGapAnSWave *fLambdaInvSq; + + ClassDef(TLambdaInvAnSWave,1) +}; + +class TLambdaInvNonMonDWave1 : public PUserFcnBase { + +public: + TLambdaInvNonMonDWave1(); + ~TLambdaInvNonMonDWave1(); + + double operator()(double, const vector&) const; + +private: + TGapNonMonDWave1 *fLambdaInvSq; + + ClassDef(TLambdaInvNonMonDWave1,1) +}; + +class TLambdaInvNonMonDWave2 : public PUserFcnBase { + +public: + TLambdaInvNonMonDWave2(); + ~TLambdaInvNonMonDWave2(); + + double operator()(double, const vector&) const; + +private: + TGapNonMonDWave2 *fLambdaInvSq; + + ClassDef(TLambdaInvNonMonDWave2,1) +}; + + +class TLambdaInvPowerLaw : public PUserFcnBase { + +public: + TLambdaInvPowerLaw() {} + ~TLambdaInvPowerLaw() {} + + double operator()(double, const vector&) const; + +private: + + ClassDef(TLambdaInvPowerLaw,1) +}; + +class TFilmMagnetizationDWave : public PUserFcnBase { + +public: + TFilmMagnetizationDWave(); + ~TFilmMagnetizationDWave(); + + double operator()(double, const vector&) const; + +private: + TGapDWave *fLambdaInvSq; + + mutable vector fPar; + + ClassDef(TFilmMagnetizationDWave,1) +}; + #endif //_TGapIntegrals_H_ diff --git a/src/external/libGapIntegrals/TGapIntegralsLinkDef.h b/src/external/libGapIntegrals/TGapIntegralsLinkDef.h index ddb102a4..005c8876 100644 --- a/src/external/libGapIntegrals/TGapIntegralsLinkDef.h +++ b/src/external/libGapIntegrals/TGapIntegralsLinkDef.h @@ -43,6 +43,19 @@ #pragma link C++ class TGapNonMonDWave1+; #pragma link C++ class TGapNonMonDWave2+; #pragma link C++ class TGapPowerLaw+; +#pragma link C++ class TLambdaSWave+; +#pragma link C++ class TLambdaDWave+; +#pragma link C++ class TLambdaAnSWave+; +#pragma link C++ class TLambdaNonMonDWave1+; +#pragma link C++ class TLambdaNonMonDWave2+; +#pragma link C++ class TLambdaPowerLaw+; +#pragma link C++ class TLambdaInvSWave+; +#pragma link C++ class TLambdaInvDWave+; +#pragma link C++ class TLambdaInvAnSWave+; +#pragma link C++ class TLambdaInvNonMonDWave1+; +#pragma link C++ class TLambdaInvNonMonDWave2+; +#pragma link C++ class TLambdaInvPowerLaw+; +#pragma link C++ class TFilmMagnetizationDWave+; #endif //__CINT__ // root dictionary stuff --------------------------------------------------