OSDN -- 오픈 소스 소프트웨어 개발 및 다운로드
Translation Status of 한국말 translation status for ko
  • R/O
  • HTTP
  • SSH
  • HTTPS
Fork

Commit

Tags
No Tags

Frequently used words (click to add to your profile)

javac++androidlinuxc#windowsobjective-ccocoa誰得qtpythonphprubygameguibathyscaphec計画中(planning stage)翻訳omegatframeworktwitterdomtestvb.netdirectxゲームエンジンbtronarduinopreviewer

Commit MetaInfo

Revision545330eb4cc38ddca93a31569d467d2d79e8b701 (tree)
Time2013-08-05 04:12:21
AuthorMikiya Fujii <mikiya.fujii@gmai...>
CommiterMikiya Fujii

Log Message

Mndo::CalcDiatomicTwoElecTwoCore1stDerivatives is refactored. #31814

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

Change Summary

Incremental Difference

--- a/src/mndo/Mndo.cpp
+++ b/src/mndo/Mndo.cpp
@@ -893,7 +893,7 @@ void Mndo::CalcCISMatrix(double** matrixCIS) const{
893893 }
894894 } // end of k-loop
895895
896- // communication to collect all matrix data on rank 0
896+ // communication to collect all matrix data on head-rank
897897 int mpiHeadRank = MolDS_mpi::MpiProcess::GetInstance()->GetHeadRank();
898898 if(mpiRank == mpiHeadRank){
899899 // receive the matrix data from other ranks
@@ -905,7 +905,7 @@ void Mndo::CalcCISMatrix(double** matrixCIS) const{
905905 }
906906 }
907907 else{
908- // send the matrix data to rank-0
908+ // send the matrix data to head-rank
909909 for(int k=0; k<this->matrixCISdimension; k++){
910910 if(k%mpiSize != mpiRank){continue;}
911911 int dest = mpiHeadRank;
@@ -969,10 +969,13 @@ double Mndo::GetCISCoefficientTwoElecIntegral(int k,
969969 return value;
970970 }
971971
972-void Mndo::MallocTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverlapAOs1stDerivs,
973- double****** diatomicOverlapAOs2ndDerivs,
974- double******* diatomicTwoElecTwoCore1stDerivs,
975- double******** diatomicTwoElecTwoCore2ndDerivs) const{
972+void Mndo::MallocTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverlapAOs1stDerivs,
973+ double****** diatomicOverlapAOs2ndDerivs,
974+ double******* diatomicTwoElecTwoCore1stDerivs,
975+ double******** diatomicTwoElecTwoCore2ndDerivs,
976+ double*** tmpRotMat,
977+ double**** tmpRotMat1stDerivs,
978+ double***** tmpDiatomicTwoElecTwoCore) const{
976979 MallocerFreer::GetInstance()->Malloc<double>(diatomicOverlapAOs1stDerivs,
977980 this->molecule->GetNumberAtoms(),
978981 OrbitalType_end,
@@ -999,12 +1002,27 @@ void Mndo::MallocTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverla
9991002 dxy,
10001003 CartesianType_end,
10011004 CartesianType_end);
1005+ MallocerFreer::GetInstance()->Malloc<double>(tmpRotMat,
1006+ OrbitalType_end,
1007+ OrbitalType_end);
1008+ MallocerFreer::GetInstance()->Malloc<double>(tmpRotMat1stDerivs,
1009+ OrbitalType_end,
1010+ OrbitalType_end,
1011+ CartesianType_end);
1012+ MallocerFreer::GetInstance()->Malloc<double>(tmpDiatomicTwoElecTwoCore,
1013+ dxy,
1014+ dxy,
1015+ dxy,
1016+ dxy);
10021017 }
10031018
10041019 void Mndo::FreeTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverlapAOs1stDerivs,
10051020 double****** diatomicOverlapAOs2ndDerivs,
10061021 double******* diatomicTwoElecTwoCore1stDerivs,
1007- double******** diatomicTwoElecTwoCore2ndDerivs) const{
1022+ double******** diatomicTwoElecTwoCore2ndDerivs,
1023+ double*** tmpRotMat,
1024+ double**** tmpRotMat1stDerivs,
1025+ double***** tmpDiatomicTwoElecTwoCore) const{
10081026 MallocerFreer::GetInstance()->Free<double>(diatomicOverlapAOs1stDerivs,
10091027 this->molecule->GetNumberAtoms(),
10101028 OrbitalType_end,
@@ -1031,6 +1049,18 @@ void Mndo::FreeTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverl
10311049 dxy,
10321050 CartesianType_end,
10331051 CartesianType_end);
1052+ MallocerFreer::GetInstance()->Free<double>(tmpRotMat,
1053+ OrbitalType_end,
1054+ OrbitalType_end);
1055+ MallocerFreer::GetInstance()->Free<double>(tmpRotMat1stDerivs,
1056+ OrbitalType_end,
1057+ OrbitalType_end,
1058+ CartesianType_end);
1059+ MallocerFreer::GetInstance()->Free<double>(tmpDiatomicTwoElecTwoCore,
1060+ dxy,
1061+ dxy,
1062+ dxy,
1063+ dxy);
10341064 }
10351065
10361066 // mu and nu is included in atomA' AO.
@@ -1572,15 +1602,21 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
15721602 stringstream ompErrors;
15731603 #pragma omp parallel
15741604 {
1575- double**** diatomicOverlapAOs1stDerivs = NULL;
1576- double***** diatomicOverlapAOs2ndDerivs = NULL;
1605+ double**** diatomicOverlapAOs1stDerivs = NULL;
1606+ double***** diatomicOverlapAOs2ndDerivs = NULL;
15771607 double****** diatomicTwoElecTwoCore1stDerivs = NULL;
15781608 double******* diatomicTwoElecTwoCore2ndDerivs = NULL;
1579- this->MallocTempMatricesEachThreadCalcHessianSCF(&diatomicOverlapAOs1stDerivs,
1580- &diatomicOverlapAOs2ndDerivs,
1581- &diatomicTwoElecTwoCore1stDerivs,
1582- &diatomicTwoElecTwoCore2ndDerivs);
1609+ double** tmpRotMat = NULL;
1610+ double*** tmpRotMat1stDerivs = NULL;
1611+ double**** tmpDiatomicTwoElecTwoCore = NULL;
15831612 try{
1613+ this->MallocTempMatricesEachThreadCalcHessianSCF(&diatomicOverlapAOs1stDerivs,
1614+ &diatomicOverlapAOs2ndDerivs,
1615+ &diatomicTwoElecTwoCore1stDerivs,
1616+ &diatomicTwoElecTwoCore2ndDerivs,
1617+ &tmpRotMat,
1618+ &tmpRotMat1stDerivs,
1619+ &tmpDiatomicTwoElecTwoCore);
15841620 #pragma omp for schedule(auto)
15851621 for(int indexAtomA=0; indexAtomA<this->molecule->GetNumberAtoms(); indexAtomA++){
15861622 const Atom& atomA = *this->molecule->GetAtom(indexAtomA);
@@ -1598,6 +1634,9 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
15981634 indexAtomA,
15991635 indexAtomB);
16001636 this->CalcDiatomicTwoElecTwoCore1stDerivatives(diatomicTwoElecTwoCore1stDerivs[indexAtomB],
1637+ tmpRotMat,
1638+ tmpRotMat1stDerivs,
1639+ tmpDiatomicTwoElecTwoCore,
16011640 indexAtomA,
16021641 indexAtomB);
16031642 this->CalcDiatomicTwoElecTwoCore2ndDerivatives(diatomicTwoElecTwoCore2ndDerivs[indexAtomB],
@@ -1659,7 +1698,10 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
16591698 this->FreeTempMatricesEachThreadCalcHessianSCF(&diatomicOverlapAOs1stDerivs,
16601699 &diatomicOverlapAOs2ndDerivs,
16611700 &diatomicTwoElecTwoCore1stDerivs,
1662- &diatomicTwoElecTwoCore2ndDerivs);
1701+ &diatomicTwoElecTwoCore2ndDerivs,
1702+ &tmpRotMat,
1703+ &tmpRotMat1stDerivs,
1704+ &tmpDiatomicTwoElecTwoCore);
16631705 }// end of omp-region
16641706 // Exception throwing for omp-region
16651707 if(!ompErrors.str().empty()){
@@ -1841,10 +1883,17 @@ void Mndo::CalcStaticFirstOrderFock(double* staticFirstOrderFock,
18411883 MallocerFreer::GetInstance()->Initialize<double>(staticFirstOrderFock,
18421884 nonRedundantQIndeces.size()+redundantQIndeces.size());
18431885 double***** diatomicTwoElecTwoCore1stDerivs = NULL;
1844- double*** diatomicOverlapAOs1stDerivs = NULL;
1886+ double*** diatomicOverlapAOs1stDerivs = NULL;
1887+ double** tmpRotMat = NULL;
1888+ double*** tmpRotMat1stDerivs = NULL;
1889+ double**** tmpDiatomicTwoElecTwoCore = NULL;
18451890
18461891 try{
1847- this->MallocTempMatricesStaticFirstOrderFock(&diatomicTwoElecTwoCore1stDerivs, &diatomicOverlapAOs1stDerivs);
1892+ this->MallocTempMatricesStaticFirstOrderFock(&diatomicTwoElecTwoCore1stDerivs,
1893+ &diatomicOverlapAOs1stDerivs,
1894+ &tmpRotMat,
1895+ &tmpRotMat1stDerivs,
1896+ &tmpDiatomicTwoElecTwoCore);
18481897 const Atom& atomA = *molecule->GetAtom(indexAtomA);
18491898 int firstAOIndexA = atomA.GetFirstAOIndex();
18501899 int lastAOIndexA = atomA.GetLastAOIndex();
@@ -1857,7 +1906,11 @@ void Mndo::CalcStaticFirstOrderFock(double* staticFirstOrderFock,
18571906 int coreChargeB = atomB.GetCoreCharge();
18581907
18591908 // calc. first derivative of two elec two core interaction
1860- this->CalcDiatomicTwoElecTwoCore1stDerivatives(diatomicTwoElecTwoCore1stDerivs, indexAtomA, indexAtomB);
1909+ this->CalcDiatomicTwoElecTwoCore1stDerivatives(diatomicTwoElecTwoCore1stDerivs,
1910+ tmpRotMat,
1911+ tmpRotMat1stDerivs,
1912+ tmpDiatomicTwoElecTwoCore,
1913+ indexAtomA, indexAtomB);
18611914 // calc. first derivative of overlapAOs.
18621915 this->CalcDiatomicOverlapAOs1stDerivatives(diatomicOverlapAOs1stDerivs, atomA, atomB);
18631916
@@ -1950,10 +2003,18 @@ void Mndo::CalcStaticFirstOrderFock(double* staticFirstOrderFock,
19502003 }
19512004 }
19522005 catch(MolDSException ex){
1953- this->FreeTempMatricesStaticFirstOrderFock(&diatomicTwoElecTwoCore1stDerivs, &diatomicOverlapAOs1stDerivs);
2006+ this->FreeTempMatricesStaticFirstOrderFock(&diatomicTwoElecTwoCore1stDerivs,
2007+ &diatomicOverlapAOs1stDerivs,
2008+ &tmpRotMat,
2009+ &tmpRotMat1stDerivs,
2010+ &tmpDiatomicTwoElecTwoCore);
19542011 throw ex;
19552012 }
1956- this->FreeTempMatricesStaticFirstOrderFock(&diatomicTwoElecTwoCore1stDerivs, &diatomicOverlapAOs1stDerivs);
2013+ this->FreeTempMatricesStaticFirstOrderFock(&diatomicTwoElecTwoCore1stDerivs,
2014+ &diatomicOverlapAOs1stDerivs,
2015+ &tmpRotMat,
2016+ &tmpRotMat1stDerivs,
2017+ &tmpDiatomicTwoElecTwoCore);
19572018
19582019 /*
19592020 printf("staticFirstOrderFock(atomA:%d axis:%s)\n",indexAtomA,CartesianTypeStr(axisA));
@@ -1964,7 +2025,10 @@ void Mndo::CalcStaticFirstOrderFock(double* staticFirstOrderFock,
19642025 }
19652026
19662027 void Mndo::MallocTempMatricesStaticFirstOrderFock(double****** diatomicTwoElecTwoCore1stDeriv,
1967- double**** diatomicOverlapAOs1stDeriv)const{
2028+ double**** diatomicOverlapAOs1stDeriv,
2029+ double*** tmpRotMat,
2030+ double**** tmpRotMat1stDerivs,
2031+ double***** tmpDiatomicTwoElecTwoCore)const{
19682032 MallocerFreer::GetInstance()->Malloc<double>(diatomicTwoElecTwoCore1stDeriv,
19692033 dxy,
19702034 dxy,
@@ -1975,10 +2039,25 @@ void Mndo::MallocTempMatricesStaticFirstOrderFock(double****** diatomicTwoElecTw
19752039 OrbitalType_end,
19762040 OrbitalType_end,
19772041 CartesianType_end);
2042+ MallocerFreer::GetInstance()->Malloc<double>(tmpRotMat,
2043+ OrbitalType_end,
2044+ OrbitalType_end);
2045+ MallocerFreer::GetInstance()->Malloc<double>(tmpRotMat1stDerivs,
2046+ OrbitalType_end,
2047+ OrbitalType_end,
2048+ CartesianType_end);
2049+ MallocerFreer::GetInstance()->Malloc<double>(tmpDiatomicTwoElecTwoCore,
2050+ dxy,
2051+ dxy,
2052+ dxy,
2053+ dxy);
19782054 }
19792055
19802056 void Mndo::FreeTempMatricesStaticFirstOrderFock(double****** diatomicTwoElecTwoCore1stDeriv,
1981- double**** diatomicOverlapAOs1stDeriv)const{
2057+ double**** diatomicOverlapAOs1stDeriv,
2058+ double*** tmpRotMat,
2059+ double**** tmpRotMat1stDerivs,
2060+ double***** tmpDiatomicTwoElecTwoCore)const{
19822061 MallocerFreer::GetInstance()->Free<double>(diatomicTwoElecTwoCore1stDeriv,
19832062 dxy,
19842063 dxy,
@@ -1989,6 +2068,18 @@ void Mndo::FreeTempMatricesStaticFirstOrderFock(double****** diatomicTwoElecTwoC
19892068 OrbitalType_end,
19902069 OrbitalType_end,
19912070 CartesianType_end);
2071+ MallocerFreer::GetInstance()->Free<double>(tmpRotMat,
2072+ OrbitalType_end,
2073+ OrbitalType_end);
2074+ MallocerFreer::GetInstance()->Free<double>(tmpRotMat1stDerivs,
2075+ OrbitalType_end,
2076+ OrbitalType_end,
2077+ CartesianType_end);
2078+ MallocerFreer::GetInstance()->Free<double>(tmpDiatomicTwoElecTwoCore,
2079+ dxy,
2080+ dxy,
2081+ dxy,
2082+ dxy);
19922083 }
19932084
19942085 // see (40) - (46) in [PT_1996].
@@ -2287,8 +2378,15 @@ void Mndo::CalcForce(const vector<int>& elecStates){
22872378 {
22882379 double***** diatomicTwoElecTwoCore1stDerivs = NULL;
22892380 double*** diatomicOverlapAOs1stDerivs = NULL;
2381+ double** tmpRotMat = NULL;
2382+ double*** tmpRotMat1stDerivs = NULL;
2383+ double**** tmpDiatomicTwoElecTwoCore = NULL;
22902384 try{
2291- this->MallocTempMatricesCalcForce(&diatomicOverlapAOs1stDerivs, &diatomicTwoElecTwoCore1stDerivs);
2385+ this->MallocTempMatricesCalcForce(&diatomicOverlapAOs1stDerivs,
2386+ &diatomicTwoElecTwoCore1stDerivs,
2387+ &tmpRotMat,
2388+ &tmpRotMat1stDerivs,
2389+ &tmpDiatomicTwoElecTwoCore);
22922390
22932391 #pragma omp for schedule(auto)
22942392 for(int b=0; b<this->molecule->GetNumberAtoms(); b++){
@@ -2300,7 +2398,11 @@ void Mndo::CalcForce(const vector<int>& elecStates){
23002398 // calc. first derivative of overlapAOs.
23012399 this->CalcDiatomicOverlapAOs1stDerivatives(diatomicOverlapAOs1stDerivs, atomA, atomB);
23022400 // calc. first derivative of two elec two core interaction
2303- this->CalcDiatomicTwoElecTwoCore1stDerivatives(diatomicTwoElecTwoCore1stDerivs, a, b);
2401+ this->CalcDiatomicTwoElecTwoCore1stDerivatives(diatomicTwoElecTwoCore1stDerivs,
2402+ tmpRotMat,
2403+ tmpRotMat1stDerivs,
2404+ tmpDiatomicTwoElecTwoCore,
2405+ a, b);
23042406
23052407 // core repulsion part
23062408 double coreRepulsion[CartesianType_end] = {0.0,0.0,0.0};
@@ -2407,7 +2509,11 @@ void Mndo::CalcForce(const vector<int>& elecStates){
24072509 #pragma omp critical
24082510 ex.Serialize(ompErrors);
24092511 }
2410- this->FreeTempMatricesCalcForce(&diatomicOverlapAOs1stDerivs, &diatomicTwoElecTwoCore1stDerivs);
2512+ this->FreeTempMatricesCalcForce(&diatomicOverlapAOs1stDerivs,
2513+ &diatomicTwoElecTwoCore1stDerivs,
2514+ &tmpRotMat,
2515+ &tmpRotMat1stDerivs,
2516+ &tmpDiatomicTwoElecTwoCore);
24112517 } // end of omp-parallelized region
24122518 // Exception throwing for omp-region
24132519 if(!ompErrors.str().empty()){
@@ -2418,10 +2524,13 @@ void Mndo::CalcForce(const vector<int>& elecStates){
24182524 // communication to reduce thsi->matrixForce on all node (namely, all_reduce)
24192525 int numTransported = elecStates.size()*this->molecule->GetNumberAtoms()*CartesianType_end;
24202526 MolDS_mpi::MpiProcess::GetInstance()->AllReduce(&this->matrixForce[0][0][0], numTransported, std::plus<double>());
2421-
24222527 }
24232528
2424-void Mndo::MallocTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs, double****** diatomicTwoElecTwoCore1stDerivs) const{
2529+void Mndo::MallocTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs,
2530+ double****** diatomicTwoElecTwoCore1stDerivs,
2531+ double*** tmpRotMat,
2532+ double**** tmpRotMat1stDerivs,
2533+ double***** tmpDiatomicTwoElecTwoCore) const{
24252534 MallocerFreer::GetInstance()->Malloc<double>(diatomicOverlapAOs1stDerivs,
24262535 OrbitalType_end,
24272536 OrbitalType_end,
@@ -2432,9 +2541,25 @@ void Mndo::MallocTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs, d
24322541 dxy,
24332542 dxy,
24342543 CartesianType_end);
2544+ MallocerFreer::GetInstance()->Malloc<double>(tmpRotMat,
2545+ OrbitalType_end,
2546+ OrbitalType_end);
2547+ MallocerFreer::GetInstance()->Malloc<double>(tmpRotMat1stDerivs,
2548+ OrbitalType_end,
2549+ OrbitalType_end,
2550+ CartesianType_end);
2551+ MallocerFreer::GetInstance()->Malloc<double>(tmpDiatomicTwoElecTwoCore,
2552+ dxy,
2553+ dxy,
2554+ dxy,
2555+ dxy);
24352556 }
24362557
2437-void Mndo::FreeTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs, double****** diatomicTwoElecTwoCore1stDerivs) const{
2558+void Mndo::FreeTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs,
2559+ double****** diatomicTwoElecTwoCore1stDerivs,
2560+ double*** tmpRotMat,
2561+ double**** tmpRotMat1stDerivs,
2562+ double***** tmpDiatomicTwoElecTwoCore) const{
24382563 MallocerFreer::GetInstance()->Free<double>(diatomicOverlapAOs1stDerivs,
24392564 OrbitalType_end,
24402565 OrbitalType_end,
@@ -2445,6 +2570,18 @@ void Mndo::FreeTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs, dou
24452570 dxy,
24462571 dxy,
24472572 CartesianType_end);
2573+ MallocerFreer::GetInstance()->Free<double>(tmpRotMat,
2574+ OrbitalType_end,
2575+ OrbitalType_end);
2576+ MallocerFreer::GetInstance()->Free<double>(tmpRotMat1stDerivs,
2577+ OrbitalType_end,
2578+ OrbitalType_end,
2579+ CartesianType_end);
2580+ MallocerFreer::GetInstance()->Free<double>(tmpDiatomicTwoElecTwoCore,
2581+ dxy,
2582+ dxy,
2583+ dxy,
2584+ dxy);
24482585 }
24492586
24502587 // see (18) in [PT_1997]
@@ -3174,6 +3311,9 @@ void Mndo::CalcDiatomicTwoElecTwoCore(double**** matrix, int indexAtomA, int ind
31743311 // Note taht d-orbital cannot be treated,
31753312 // that is, matrix[dxy][dxy][dxy][dxy][CartesianType_end] cannot be treatable.
31763313 void Mndo::CalcDiatomicTwoElecTwoCore1stDerivatives(double***** matrix,
3314+ double** tmpRotMat,
3315+ double*** tmpRotMat1stDerivs,
3316+ double**** tmpDiatomicTwoElecTwoCore,
31773317 int indexAtomA,
31783318 int indexAtomB) const{
31793319 const Atom& atomA = *this->molecule->GetAtom(indexAtomA);
@@ -3200,65 +3340,48 @@ void Mndo::CalcDiatomicTwoElecTwoCore1stDerivatives(double***** matrix,
32003340 dxy,
32013341 CartesianType_end);
32023342
3203- double** rotatingMatrix = NULL;
3204- double*** rotMat1stDerivatives = NULL;
3205- double**** diatomicTwoElecTwoCore = NULL;
3206- try{
3207- this->MallocDiatomicTwoElecTwoCore1stDeriTemps(&rotatingMatrix,
3208- &rotMat1stDerivatives,
3209- &diatomicTwoElecTwoCore);
3210- // calclation in diatomic frame
3211- for(int mu=0; mu<atomA.GetValenceSize(); mu++){
3212- for(int nu=mu; nu<atomA.GetValenceSize(); nu++){
3213- for(int lambda=0; lambda<atomB.GetValenceSize(); lambda++){
3214- for(int sigma=lambda; sigma<atomB.GetValenceSize(); sigma++){
3215- for(int dimA=0; dimA<CartesianType_end; dimA++){
3216- matrix[mu][nu][lambda][sigma][dimA]
3217- = this->GetNddoRepulsionIntegral1stDerivative(
3218- atomA,
3219- atomA.GetValence(mu),
3220- atomA.GetValence(nu),
3221- atomB,
3222- atomB.GetValence(lambda),
3223- atomB.GetValence(sigma),
3224- static_cast<CartesianType>(dimA));
3225- matrix[nu][mu][lambda][sigma][dimA] = matrix[mu][nu][lambda][sigma][dimA];
3226- matrix[nu][mu][sigma][lambda][dimA] = matrix[mu][nu][lambda][sigma][dimA];
3227- matrix[mu][nu][sigma][lambda][dimA] = matrix[mu][nu][lambda][sigma][dimA];
3228- }
3229- diatomicTwoElecTwoCore[mu][nu][lambda][sigma]
3230- = this->GetNddoRepulsionIntegral(
3343+ // calclation in diatomic frame
3344+ for(int mu=0; mu<atomA.GetValenceSize(); mu++){
3345+ for(int nu=mu; nu<atomA.GetValenceSize(); nu++){
3346+ for(int lambda=0; lambda<atomB.GetValenceSize(); lambda++){
3347+ for(int sigma=lambda; sigma<atomB.GetValenceSize(); sigma++){
3348+ for(int dimA=0; dimA<CartesianType_end; dimA++){
3349+ matrix[mu][nu][lambda][sigma][dimA]
3350+ = this->GetNddoRepulsionIntegral1stDerivative(
32313351 atomA,
32323352 atomA.GetValence(mu),
32333353 atomA.GetValence(nu),
32343354 atomB,
32353355 atomB.GetValence(lambda),
3236- atomB.GetValence(sigma));
3237- diatomicTwoElecTwoCore[nu][mu][lambda][sigma] = diatomicTwoElecTwoCore[mu][nu][lambda][sigma];
3238- diatomicTwoElecTwoCore[nu][mu][sigma][lambda] = diatomicTwoElecTwoCore[mu][nu][lambda][sigma];
3239- diatomicTwoElecTwoCore[mu][nu][sigma][lambda] = diatomicTwoElecTwoCore[mu][nu][lambda][sigma];
3240- }
3356+ atomB.GetValence(sigma),
3357+ static_cast<CartesianType>(dimA));
3358+ matrix[nu][mu][lambda][sigma][dimA] = matrix[mu][nu][lambda][sigma][dimA];
3359+ matrix[nu][mu][sigma][lambda][dimA] = matrix[mu][nu][lambda][sigma][dimA];
3360+ matrix[mu][nu][sigma][lambda][dimA] = matrix[mu][nu][lambda][sigma][dimA];
3361+ }
3362+ tmpDiatomicTwoElecTwoCore[mu][nu][lambda][sigma]
3363+ = this->GetNddoRepulsionIntegral(
3364+ atomA,
3365+ atomA.GetValence(mu),
3366+ atomA.GetValence(nu),
3367+ atomB,
3368+ atomB.GetValence(lambda),
3369+ atomB.GetValence(sigma));
3370+ tmpDiatomicTwoElecTwoCore[nu][mu][lambda][sigma] = tmpDiatomicTwoElecTwoCore[mu][nu][lambda][sigma];
3371+ tmpDiatomicTwoElecTwoCore[nu][mu][sigma][lambda] = tmpDiatomicTwoElecTwoCore[mu][nu][lambda][sigma];
3372+ tmpDiatomicTwoElecTwoCore[mu][nu][sigma][lambda] = tmpDiatomicTwoElecTwoCore[mu][nu][lambda][sigma];
32413373 }
32423374 }
32433375 }
3244-
3245- // rotate matirix into the space frame
3246- this->CalcRotatingMatrix(rotatingMatrix, atomA, atomB);
3247- this->CalcRotatingMatrix1stDerivatives(rotMat1stDerivatives, atomA, atomB);
3248- this->RotateDiatomicTwoElecTwoCore1stDerivativesToSpaceFrame(matrix,
3249- diatomicTwoElecTwoCore,
3250- rotatingMatrix,
3251- rotMat1stDerivatives);
32523376 }
3253- catch(MolDSException ex){
3254- this->FreeDiatomicTwoElecTwoCore1stDeriTemps(&rotatingMatrix,
3255- &rotMat1stDerivatives,
3256- &diatomicTwoElecTwoCore);
3257- throw ex;
3258- }
3259- this->FreeDiatomicTwoElecTwoCore1stDeriTemps(&rotatingMatrix,
3260- &rotMat1stDerivatives,
3261- &diatomicTwoElecTwoCore);
3377+
3378+ // rotate matirix into the space frame
3379+ this->CalcRotatingMatrix(tmpRotMat, atomA, atomB);
3380+ this->CalcRotatingMatrix1stDerivatives(tmpRotMat1stDerivs, atomA, atomB);
3381+ this->RotateDiatomicTwoElecTwoCore1stDerivativesToSpaceFrame(matrix,
3382+ tmpDiatomicTwoElecTwoCore,
3383+ tmpRotMat,
3384+ tmpRotMat1stDerivs);
32623385 }
32633386
32643387 // Calculation of second derivatives of the two electrons two cores integral in space fixed frame,
--- a/src/mndo/Mndo.h
+++ b/src/mndo/Mndo.h
@@ -150,11 +150,17 @@ private:
150150 void MallocTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverlapAOs1stDerivs,
151151 double****** diatomicOverlapAOs2ndDerivs,
152152 double******* diatomicTwoElecTwoCore1stDerivs,
153- double******** diatomicTwoElecTwoCore2ndDerivs) const;
153+ double******** diatomicTwoElecTwoCore2ndDerivs,
154+ double*** tmpRotMat,
155+ double**** tmpRotMat1stDerivs,
156+ double***** tmpDiatomicTwoElecTwo) const;
154157 void FreeTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverlapAOs1stDerivs,
155158 double****** diatomicOverlapAOs2ndDerivs,
156- double******* diatomicTwoElecTwoCore1stDeriv,
157- double******** diatomicTwoElecTwoCore2ndDerivs) const;
159+ double******* diatomicTwoElecTwoCore1stDerivs,
160+ double******** diatomicTwoElecTwoCore2ndDerivs,
161+ double*** tmpRotMat,
162+ double**** tmpRotMat1stDerivs,
163+ double***** tmpDiatomicTwoElecTwo) const;
158164 double GetAuxiliaryHessianElement1(int mu,
159165 int nu,
160166 int indexAtomA,
@@ -241,9 +247,15 @@ private:
241247 int indexAtomA,
242248 MolDS_base::CartesianType axisA) const;
243249 void MallocTempMatricesStaticFirstOrderFock(double****** diatomicTwoElecTwoCore1stDeriv,
244- double**** diatomicOverlapAOs1stDeriv) const;
250+ double**** diatomicOverlapAOs1stDeriv,
251+ double*** tmpRotMat,
252+ double**** tmpRotMat1stDerivs,
253+ double***** tmpDiatomicTwoElecTwo) const;
245254 void FreeTempMatricesStaticFirstOrderFock(double****** diatomicTwoElecTwoCore1stDeriv,
246- double**** diatomicOverlapAOs1stDeriv) const;
255+ double**** diatomicOverlapAOs1stDeriv,
256+ double*** tmpRotMat,
257+ double**** tmpRotMat1stDerivs,
258+ double***** tmpDiatomicTwoElecTwo) const;
247259 void CalcMatrixCPHF(double** matrixCPHF,
248260 const std::vector<MoIndexPair>& nonRedundantQIndeces,
249261 const std::vector<MoIndexPair>& redundantQIndeces) const;
@@ -266,6 +278,9 @@ private:
266278 MolDS_base::CartesianType axisA) const;
267279 void CalcDiatomicTwoElecTwoCore(double**** matrix, int indexAtomA, int indexAtomB) const;
268280 void CalcDiatomicTwoElecTwoCore1stDerivatives(double***** matrix,
281+ double** tmpRotMat,
282+ double*** tmpRotMat1stDerivs,
283+ double**** tmpDiatomicTwoElecTwoCore,
269284 int indexAtomA,
270285 int indexAtomB) const;
271286 void CalcDiatomicTwoElecTwoCore2ndDerivatives(double****** matrix,
@@ -356,9 +371,15 @@ private:
356371 MolDS_base::MultipoleType multipoleB,
357372 double rAB) const;
358373 void MallocTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs,
359- double****** diatomiTwoElecTwoCore1stDerivs) const;
374+ double****** diatomiTwoElecTwoCore1stDerivs,
375+ double*** tmpRotMat,
376+ double**** tmpRotMat1stDerivs,
377+ double***** tmpDiatomicTwoElecTwoCore) const;
360378 void FreeTempMatricesCalcForce(double**** diatomicOverlapAOs1stDerivs,
361- double****** diatomiTwoElecTwoCore1stDerivs) const;
379+ double****** diatomiTwoElecTwoCore1stDerivs,
380+ double*** tmpRotMat,
381+ double**** tmpRotMat1stDerivs,
382+ double***** tmpDiatomicTwoElecTwoCore) const;
362383 void CalcForceSCFElecCoreAttractionPart(double* force,
363384 int indexAtomA,
364385 int indexAtomB,