Commit: cd8563de69c739f1a7371700b6fce69091e3e18c - MolDS (git)

Commit MetaInfo

• Revisión: cd8563de69c739f1a7371700b6fce69091e3e18c (tree)
• Tiempo: 2013-01-22 00:11:44
• Autor: okymmolds <okymmolds@gmai...>
• Commiter: okymmolds

Log Message

Unpaired electron population at Mulliken type. #29958

git-svn-id: https://svn.sourceforge.jp/svnroot/molds/trunk@1236 1136aad2-a195-0410-b898-f5ea1d11b9d8

Cambiar Resumen

• modified: src/base/InputParser.cpp
• modified: src/base/InputParser.h
• modified: src/base/Parameters.cpp
• modified: src/base/Parameters.h
• modified: src/cndo/Cndo2.cpp
• modified: src/cndo/Cndo2.h
• modified: src/zindo/ZindoS.cpp
• modified: src/zindo/ZindoS.h

Diferencia incremental

--- a/src/base/InputParser.cpp
+++ b/src/base/InputParser.cpp
@@ -311,6 +311,7 @@ void InputParser::SetMessages(){
    this->stringCISAllTransitionDipoleMoments = "all_transition_dipole_moments";
    this->stringCISNumPrintCoefficients       = "num_print_coefficients";
    this->stringCISMulliken                   = "mulliken";
+   this->stringCISUnpairedPop                = "unpaired_electron_population";
 
    // Memory
    this->stringMemory          = "memory";
@@ -823,6 +824,16 @@ int InputParser::ParseConditionsCIS(vector<string>* inputTerms, int parseIndex)
          }
          parseIndex++;
       }
+      // unpaired electron population
+      if((*inputTerms)[parseIndex].compare(this->stringCISUnpairedPop) == 0){
+          if((*inputTerms)[parseIndex+1].compare(this->stringYES) == 0){
+             Parameters::GetInstance()->SetRequiresUnpairedPopCIS(true);
+          }
+          else {
+             Parameters::GetInstance()->SetRequiresUnpairedPopCIS(false);
+          }
+          parseIndex++;   
+      }
       parseIndex++;   
    }
    return parseIndex;

--- a/src/base/InputParser.h
+++ b/src/base/InputParser.h
@@ -236,6 +236,7 @@ private:
    std::string stringCISAllTransitionDipoleMoments;
    std::string stringCISNumPrintCoefficients;
    std::string stringCISMulliken;
+   std::string stringCISUnpairedPop;
    // Memory
    std::string stringMemory;
    std::string stringMemoryEnd;

--- a/src/base/Parameters.cpp
+++ b/src/base/Parameters.cpp
@@ -1,6 +1,7 @@
 //************************************************************************//
 // Copyright (C) 2011-2012 Mikiya Fujii                                   // 
 // Copyright (C) 2012-2012 Katsuhiko Nishimra                             // 
+// Copyright (C) 2012-2013 Michihiro Okuyama
 //                                                                        // 
 // This file is part of MolDS.                                            // 
 //                                                                        // 
@@ -171,6 +172,7 @@ void Parameters::SetDefaultValues(){
    this->numberPrintCoefficientsCIS            = 1;
    this->requiresExcitonEnergiesCIS            = false;
    this->requiresAllTransitionDipoleMomentsCIS = false;
+   this->requiresUnpairedPopCIS                = false;
    // Memory
    this->limitHeapMemory = 256;
    // MD
@@ -705,6 +707,14 @@ bool Parameters::RequiresMullikenCIS() const{
            0<this->electronicStateIndecesMullikenCIS->size());
 }
 
+void Parameters::SetRequiresUnpairedPopCIS(bool requiresUnpairedPopCIS){
+   this->requiresUnpairedPopCIS = requiresUnpairedPopCIS;
+}
+
+bool Parameters::RequiresUnpairedPopCIS() const{
+   return this->requiresUnpairedPopCIS; 
+}
+
 // methods for Memory
 double Parameters::GetLimitHeapMemory() const{
    return this->limitHeapMemory;

--- a/src/base/Parameters.h
+++ b/src/base/Parameters.h
@@ -1,6 +1,7 @@
 //************************************************************************//
 // Copyright (C) 2011-2012 Mikiya Fujii                                   // 
 // Copyright (C) 2012-2012 Katsuhiko Nishimra                             // 
+// Copyright (C) 2012-2013 Michihiro Okuyama
 //                                                                        // 
 // This file is part of MolDS.                                            // 
 //                                                                        // 
@@ -141,6 +142,8 @@ public:
    std::vector<int>* GetElectronicStateIndecesMullikenCIS() const;
    void              AddElectronicStateIndexMullikenCIS(int electronicStateIndex);
    bool              RequiresMullikenCIS() const;
+   void              SetRequiresUnpairedPopCIS(bool requiresUnpairedPopCIS); 
+   bool              RequiresUnpairedPopCIS() const; 
    // Memory
    double GetLimitHeapMemory() const;
    void   SetLimitHeapMemory(double limitHeap);
@@ -285,6 +288,7 @@ private:
    bool              requiresExcitonEnergiesCIS;
    bool              requiresAllTransitionDipoleMomentsCIS;
    std::vector<int>* electronicStateIndecesMullikenCIS;
+   bool              requiresUnpairedPopCIS;
    // Memory
    double limitHeapMemory;
    // MD

--- a/src/cndo/Cndo2.cpp
+++ b/src/cndo/Cndo2.cpp
@@ -1,5 +1,6 @@
 //************************************************************************//
 // Copyright (C) 2011-2012 Mikiya Fujii                                   // 
+// Copyright (C) 2012-2013 Michihiro Okuyama
 //                                                                        // 
 // This file is part of MolDS.                                            // 
 //                                                                        // 
@@ -72,6 +73,7 @@ Cndo2::Cndo2(){
    this->orbitalElectronPopulationCIS = NULL;
    this->atomicElectronPopulation     = NULL;
    this->atomicElectronPopulationCIS  = NULL;
+   this->atomicUnpairedPopulationCIS  = NULL;
    this->overlapAOs = NULL;
    this->twoElecTwoCore = NULL;
    this->cartesianMatrix = NULL;
@@ -205,6 +207,8 @@ void Cndo2::SetMessages(){
    this->messageMullikenAtomsSCF   = "\tMulliken charge(SCF):";
    this->messageMullikenAtoms      = "\tMulliken charge:";
    this->messageMullikenAtomsTitle = "\t\t\t\t| k-th eigenstate | i-th atom | atom type | core charge[a.u.] | Mulliken charge[a.u.]| \n";
+   this->messageUnpairedAtoms      = "\tUnpaired electron population:";
+   this->messageUnpairedAtomsTitle = "\t\t\t\t| k-th eigenstate | i-th atom | atom type | Unpaired electron population[a.u.]| \n";
    this->messageElecEnergy = "\tElectronic energy(SCF):";
    this->messageNoteElecEnergy = "\tNote that this electronic energy includes core-repulsions.\n\n";
    this->messageNoteElecEnergyVdW = "\tNote that this electronic energy includes core-repulsions and vdW correction.\n\n";

--- a/src/cndo/Cndo2.h
+++ b/src/cndo/Cndo2.h
@@ -1,5 +1,6 @@
 //************************************************************************//
 // Copyright (C) 2011-2012 Mikiya Fujii                                   // 
+// Copyright (C) 2012-2013 Michihiro Okuyama
 //                                                                        // 
 // This file is part of MolDS.                                            // 
 //                                                                        // 
@@ -86,6 +87,8 @@ protected:
    std::string messageOmpElapsedTimeSCF;
    std::string messageMullikenAtoms;
    std::string messageMullikenAtomsTitle;
+   std::string messageUnpairedAtoms;
+   std::string messageUnpairedAtomsTitle;
    std::string messageUnitSec; 
    std::vector<MolDS_base::AtomType> enableAtomTypes;
    MolDS_base::Molecule* molecule;
@@ -99,6 +102,7 @@ protected:
    double***    orbitalElectronPopulationCIS; 
    double*      atomicElectronPopulation; //P_{AB} of (3.21) in J. A. Pople book.
    double**     atomicElectronPopulationCIS; 
+   double**     atomicUnpairedPopulationCIS; 
    double**     overlapAOs; // overlap integral between AOs
    double****** twoElecTwoCore;
    double***    cartesianMatrix; // cartesian matrix represented by AOs

--- a/src/zindo/ZindoS.cpp
+++ b/src/zindo/ZindoS.cpp
@@ -1,5 +1,6 @@
 //************************************************************************//
 // Copyright (C) 2011-2012 Mikiya Fujii                                   // 
+// Copyright (C) 2012-2013 Michihiro Okuyama
 //                                                                        // 
 // This file is part of MolDS.                                            // 
 //                                                                        // 
@@ -1020,6 +1021,10 @@ void ZindoS::CalcCISProperties(){
    this->CalcAtomicElectronPopulationCIS(&this->atomicElectronPopulationCIS,
                                          this->orbitalElectronPopulationCIS, 
                                          *this->molecule);
+   // atomic unpaired electron population
+   this->CalcAtomicUnpairedPopulationCIS(&this->atomicUnpairedPopulationCIS,
+                                         this->orbitalElectronPopulationCIS, 
+                                         *this->molecule);
 }
 
 void ZindoS::CalcElectronicDipoleMomentsExcitedState(double*** electronicTransitionDipoleMoments,
@@ -1406,6 +1411,58 @@ void ZindoS::CalcAtomicElectronPopulationCIS(double*** atomicElectronPopulationC
    }
 }
 
+void ZindoS::CalcAtomicUnpairedPopulationCIS(double*** atomicUnpairedPopulationCIS,
+                                             double const* const* const* orbitalElectronPopulationCIS, 
+                                             const Molecule& molecule) const{
+   if(!Parameters::GetInstance()->RequiresMullikenCIS()){
+      return;
+   }
+   if(!Parameters::GetInstance()->RequiresUnpairedPopCIS()){
+      return;
+   }
+   int totalNumberAtoms = molecule.GetNumberAtoms();
+   vector<int>* elecStates = Parameters::GetInstance()->GetElectronicStateIndecesMullikenCIS();
+   // malloc or initialize free exciton energies
+   if(*atomicUnpairedPopulationCIS == NULL){
+      MallocerFreer::GetInstance()->Malloc<double>(atomicUnpairedPopulationCIS, 
+                                                   elecStates->size(),
+                                                   totalNumberAtoms);
+   }
+   else{
+      MallocerFreer::GetInstance()->Initialize<double>(*atomicUnpairedPopulationCIS,
+                                                       elecStates->size(),
+                                                       totalNumberAtoms);
+   }
+   // calc atomic electron population
+   for(int k=0; k<elecStates->size(); k++){
+      stringstream ompErrors;
+#pragma omp parallel for schedule(auto)
+      for(int a=0; a<totalNumberAtoms; a++){
+         try{
+            int firstAOIndex = molecule.GetAtom(a)->GetFirstAOIndex();
+            int numberAOs = molecule.GetAtom(a)->GetValenceSize();
+            (*atomicUnpairedPopulationCIS)[k][a] = 0.0; 
+            for(int i=firstAOIndex; i<firstAOIndex+numberAOs; i++){
+               double orbitalSquarePopulation = 0.0; 
+               int    totalNumberAOs = molecule.GetTotalNumberAOs(); 
+               for(int j=0; j<totalNumberAOs; j++) {
+                  orbitalSquarePopulation += orbitalElectronPopulationCIS[k][i][j] * orbitalElectronPopulationCIS[k][j][i]; 
+               }
+               (*atomicUnpairedPopulationCIS)[k][a] += 2.0 * orbitalElectronPopulationCIS[k][i][i] - orbitalSquarePopulation;
+            }
+         }
+         catch(MolDSException ex){
+#pragma omp critical
+            ompErrors << ex.what() << endl ;
+         }
+      }
+      // Exception throwing for omp-region
+      if(!ompErrors.str().empty()){
+         throw MolDSException(ompErrors.str());
+      }
+   }
+}
+
 void ZindoS::OutputCISResults() const{
    int numberActiveOcc = Parameters::GetInstance()->GetActiveOccCIS();
    int numberActiveVir = Parameters::GetInstance()->GetActiveVirCIS();
@@ -1441,6 +1498,9 @@ void ZindoS::OutputCISResults() const{
    // output mulliken population
    this->OutputCISMulliken();
 
+   // output unpaired electron population
+   this->OutputCISUnpairedPop(); 
+
    // output exciton energies
    if(Parameters::GetInstance()->RequiresExcitonEnergiesCIS()){
       this->OutputLog(this->messageExcitonEnergiesCIS);
@@ -1588,6 +1648,29 @@ void ZindoS::OutputCISMulliken() const{
    }
 }
 
+void ZindoS::OutputCISUnpairedPop() const{
+   if(!Parameters::GetInstance()->RequiresMullikenCIS()){
+      return;
+   }
+   if(!Parameters::GetInstance()->RequiresUnpairedPopCIS()){
+      return;
+   }
+   int totalNumberAtoms = this->molecule->GetNumberAtoms();
+   this->OutputLog(this->messageUnpairedAtomsTitle);
+   vector<int>* elecStates = Parameters::GetInstance()->GetElectronicStateIndecesMullikenCIS();
+   for(int k=0; k<elecStates->size(); k++){
+      for(int a=0; a<totalNumberAtoms; a++){
+         const Atom& atom = *this->molecule->GetAtom(a);
+         this->OutputLog(boost::format("%s\t%d\t%d\t%s\t%e\n") % this->messageUnpairedAtoms
+                                                               % (*elecStates)[k]
+                                                               % a
+                                                               % AtomTypeStr(atom.GetAtomType())
+                                                               % this->atomicUnpairedPopulationCIS[k][a]);
+      }
+      this->OutputLog("\n");
+   }
+}
+
 void ZindoS::SortCISEigenVectorCoefficients(vector<CISEigenVectorCoefficient>* cisEigenVectorCoefficients,
                                             double* cisEigenVector) const{
    for(int l=0; l<this->matrixCISdimension; l++){

--- a/src/zindo/ZindoS.h
+++ b/src/zindo/ZindoS.h
@@ -1,5 +1,6 @@
 //************************************************************************//
 // Copyright (C) 2011-2012 Mikiya Fujii                                   // 
+// Copyright (C) 2012-2013 Michihiro Okuyama
 //                                                                        // 
 // This file is part of MolDS.                                            // 
 //                                                                        // 
@@ -139,6 +140,7 @@ private:
    void OutputCISDipole() const;
    void OutputCISTransitionDipole() const;
    void OutputCISMulliken() const;
+   void OutputCISUnpairedPop() const;
    void CalcFreeExcitonEnergies(double** freeExcitonEnergiesCIS, 
                                 const MolDS_base::Molecule& molecule, 
                                 double const* energiesMO, 
@@ -152,6 +154,9 @@ private:
    void CalcAtomicElectronPopulationCIS(double*** atomicElectronPopulationCIS,
                                         double const* const* const* orbitalElectronPopulationCIS, 
                                         const MolDS_base::Molecule& molecule) const;
+   void CalcAtomicUnpairedPopulationCIS(double*** atomicUnpairedPopulationCIS,
+                                        double const* const* const* orbitalElectronPopulationCIS, 
+                                        const MolDS_base::Molecule& molecule) const; 
    void CalcElectronicDipoleMomentsExcitedState(double*** electronicTransitionDipoleMoments,
                                                 double const* const* fockMatrix,
                                                 double const* const* matrixCIS,