• R/O
  • HTTP
  • SSH
  • HTTPS

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

Revisãoe96cf5382f3057d6f4991e58ffe16f3fbc4a182d (tree)
Hora2013-08-06 18:16:21
AutorMikiya Fujii <mikiya.fujii@gmai...>
CommiterMikiya Fujii

Mensagem de Log

Cndo2::CalcDiatomicOverlapAOs2ndDerivatives is refactored also. #31824

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

Mudança Sumário

Diff

--- a/src/cndo/Cndo2.cpp
+++ b/src/cndo/Cndo2.cpp
@@ -3912,6 +3912,14 @@ void Cndo2::CalcDiatomicOverlapAOs1stDerivatives(double*** diatomicOverlapAOs1st
39123912 // Note that this method can not treat d-obitals
39133913 // because CalcRotatingMatrix1stDerivatives can not treat d-orbitals.
39143914 void Cndo2::CalcDiatomicOverlapAOs2ndDerivatives(double**** diatomicOverlapAOs2ndDerivs,
3915+ double** tmpDiaOverlapAOsInDiaFrame,
3916+ double** tmpDiaOverlapAOs1stDerivInDiaFrame,
3917+ double** tmpDiaOverlapAOs2ndDerivInDiaFrame,
3918+ double*** tmpDiaOverlapAOs1stDerivs,
3919+ double**** tmpDiaOverlapAOs2ndDerivs,
3920+ double** tmpRotMat,
3921+ double*** tmpRotMat1stDerivs,
3922+ double**** tmpRotMat2ndDerivs,
39153923 const Atom& atomA,
39163924 const Atom& atomB) const{
39173925 double cartesian[CartesianType_end] = {atomA.GetXyz()[XAxis] - atomB.GetXyz()[XAxis],
@@ -3921,129 +3929,90 @@ void Cndo2::CalcDiatomicOverlapAOs2ndDerivatives(double**** diatomicOverlapAOs2n
39213929 pow(cartesian[YAxis],2.0) +
39223930 pow(cartesian[ZAxis],2.0) );
39233931
3924- double** diaOverlapAOsInDiaFrame = NULL; // diatomic overlapAOs in diatomic frame
3925- double** diaOverlapAOs1stDerivInDiaFrame = NULL; // first derivative of the diaOverlapAOs. This derivative is related to the distance between two atoms.
3926- double** diaOverlapAOs2ndDerivInDiaFrame = NULL; // second derivative of the diaOverlapAOs. This derivative is related to the distance between two atoms.
3927- double** rotMat = NULL; // rotating Matrix from the diatomic frame to space fixed frame.
3928- double*** rotMat1stDerivatives = NULL; //first derivatives of the rotMat
3929- double**** rotMat2ndDerivatives = NULL; //second derivatives of the rotMat
3930- double*** tempDiaOverlapAOs1stDerivs = NULL; // first derivatives of the diaOverlapAOs. This derivatives are related to the all Cartesian coordinates.
3931- double**** tempDiaOverlapAOs2ndDerivs = NULL; //sedond derivatives of the diaOverlapAOs. This derivatives are related to the all Cartesian coordinates.
3932+ this->CalcDiatomicOverlapAOsInDiatomicFrame(tmpDiaOverlapAOsInDiaFrame, atomA, atomB);
3933+ this->CalcDiatomicOverlapAOs1stDerivativeInDiatomicFrame(tmpDiaOverlapAOs1stDerivInDiaFrame, atomA, atomB);
3934+ this->CalcDiatomicOverlapAOs2ndDerivativeInDiatomicFrame(tmpDiaOverlapAOs2ndDerivInDiaFrame, atomA, atomB);
3935+ this->CalcRotatingMatrix(tmpRotMat, atomA, atomB);
3936+ this->CalcRotatingMatrix1stDerivatives(tmpRotMat1stDerivs, atomA, atomB);
3937+ this->CalcRotatingMatrix2ndDerivatives(tmpRotMat2ndDerivs, atomA, atomB);
39323938
3933- try{
3934- this->MallocDiatomicOverlapAOs2ndDeriTemps(&diaOverlapAOsInDiaFrame,
3935- &diaOverlapAOs1stDerivInDiaFrame,
3936- &diaOverlapAOs2ndDerivInDiaFrame,
3937- &rotMat,
3938- &rotMat1stDerivatives,
3939- &rotMat2ndDerivatives,
3940- &tempDiaOverlapAOs1stDerivs,
3941- &tempDiaOverlapAOs2ndDerivs);
3942- this->CalcDiatomicOverlapAOsInDiatomicFrame(diaOverlapAOsInDiaFrame, atomA, atomB);
3943- this->CalcDiatomicOverlapAOs1stDerivativeInDiatomicFrame(diaOverlapAOs1stDerivInDiaFrame, atomA, atomB);
3944- this->CalcDiatomicOverlapAOs2ndDerivativeInDiatomicFrame(diaOverlapAOs2ndDerivInDiaFrame, atomA, atomB);
3945- this->CalcRotatingMatrix(rotMat, atomA, atomB);
3946- this->CalcRotatingMatrix1stDerivatives(rotMat1stDerivatives, atomA, atomB);
3947- this->CalcRotatingMatrix2ndDerivatives(rotMat2ndDerivatives, atomA, atomB);
3948-
3949- // calculate each element of first derivatives
3950- for(int i=0; i<OrbitalType_end; i++){
3951- for(int j=0; j<OrbitalType_end; j++){
3952- for(int dimA1=0; dimA1<CartesianType_end; dimA1++){
3953- tempDiaOverlapAOs1stDerivs[i][j][dimA1] = (cartesian[dimA1]/R)*diaOverlapAOs1stDerivInDiaFrame[i][j];
3954- }
3939+ // calculate each element of first derivatives
3940+ for(int i=0; i<OrbitalType_end; i++){
3941+ for(int j=0; j<OrbitalType_end; j++){
3942+ for(int dimA1=0; dimA1<CartesianType_end; dimA1++){
3943+ tmpDiaOverlapAOs1stDerivs[i][j][dimA1] = (cartesian[dimA1]/R)*tmpDiaOverlapAOs1stDerivInDiaFrame[i][j];
39553944 }
39563945 }
3946+ }
39573947
3958- // calculate each element of second derivatives
3959- for(int i=0; i<OrbitalType_end; i++){
3960- for(int j=0; j<OrbitalType_end; j++){
3961- for(int dimA1=XAxis; dimA1<CartesianType_end; dimA1++){
3962- for(int dimA2=XAxis; dimA2<CartesianType_end; dimA2++){
3963- tempDiaOverlapAOs2ndDerivs[i][j][dimA1][dimA2]
3964- = this->Get2ndDerivativeElementFromDistanceDerivatives(diaOverlapAOs1stDerivInDiaFrame[i][j],
3965- diaOverlapAOs2ndDerivInDiaFrame[i][j],
3966- static_cast<CartesianType>(dimA1),
3967- static_cast<CartesianType>(dimA2),
3968- cartesian,
3969- R);
3970- }
3948+ // calculate each element of second derivatives
3949+ for(int i=0; i<OrbitalType_end; i++){
3950+ for(int j=0; j<OrbitalType_end; j++){
3951+ for(int dimA1=XAxis; dimA1<CartesianType_end; dimA1++){
3952+ for(int dimA2=XAxis; dimA2<CartesianType_end; dimA2++){
3953+ tmpDiaOverlapAOs2ndDerivs[i][j][dimA1][dimA2]
3954+ = this->Get2ndDerivativeElementFromDistanceDerivatives(tmpDiaOverlapAOs1stDerivInDiaFrame[i][j],
3955+ tmpDiaOverlapAOs2ndDerivInDiaFrame[i][j],
3956+ static_cast<CartesianType>(dimA1),
3957+ static_cast<CartesianType>(dimA2),
3958+ cartesian,
3959+ R);
39713960 }
39723961 }
39733962 }
3963+ }
39743964
3975- // rotate
3976- for(int i=0; i<OrbitalType_end; i++){
3977- for(int j=0; j<OrbitalType_end; j++){
3978- for(int dimA1=XAxis; dimA1<CartesianType_end; dimA1++){
3979- for(int dimA2=XAxis; dimA2<CartesianType_end; dimA2++){
3980- diatomicOverlapAOs2ndDerivs[i][j][dimA1][dimA2] = 0.0;
3981-
3982- double temp1 = 0.0, temp2=0.0, temp3 = 0.0;
3983- double temp4 = 0.0, temp5=0.0, temp6 = 0.0;
3984- double temp7 = 0.0, temp8=0.0, temp9 = 0.0;
3985- for(int k=0; k<OrbitalType_end; k++){
3986- for(int l=0; l<OrbitalType_end; l++){
3987-
3988- temp1 += rotMat2ndDerivatives [i][k][dimA1][dimA2]
3989- *rotMat [j][l]
3990- *diaOverlapAOsInDiaFrame[k][l];
3991- temp2 += rotMat [i][k]
3992- *rotMat2ndDerivatives [j][l][dimA1][dimA2]
3993- *diaOverlapAOsInDiaFrame[k][l];
3994- temp3 += rotMat [i][k]
3995- *rotMat [j][l]
3996- *tempDiaOverlapAOs2ndDerivs[k][l][dimA1][dimA2];
3997- temp4 += rotMat1stDerivatives [i][k][dimA1]
3998- *rotMat1stDerivatives [j][l][dimA2]
3999- *diaOverlapAOsInDiaFrame[k][l];
4000- temp5 += rotMat1stDerivatives [i][k][dimA1]
4001- *rotMat [j][l]
4002- *tempDiaOverlapAOs1stDerivs[k][l][dimA2];
4003- temp6 += rotMat1stDerivatives [i][k][dimA2]
4004- *rotMat1stDerivatives [j][l][dimA1]
4005- *diaOverlapAOsInDiaFrame[k][l];
4006- temp7 += rotMat [i][k]
4007- *rotMat1stDerivatives [j][l][dimA1]
4008- *tempDiaOverlapAOs1stDerivs[k][l][dimA2];
4009- temp8 += rotMat1stDerivatives [i][k][dimA2]
4010- *rotMat [j][l]
4011- *tempDiaOverlapAOs1stDerivs[k][l][dimA1];
4012- temp9 += rotMat [i][k]
4013- *rotMat1stDerivatives [j][l][dimA2]
4014- *tempDiaOverlapAOs1stDerivs[k][l][dimA1];
4015- }
3965+ // rotate
3966+ for(int i=0; i<OrbitalType_end; i++){
3967+ for(int j=0; j<OrbitalType_end; j++){
3968+ for(int dimA1=XAxis; dimA1<CartesianType_end; dimA1++){
3969+ for(int dimA2=XAxis; dimA2<CartesianType_end; dimA2++){
3970+ diatomicOverlapAOs2ndDerivs[i][j][dimA1][dimA2] = 0.0;
3971+
3972+ double temp1 = 0.0, temp2=0.0, temp3 = 0.0;
3973+ double temp4 = 0.0, temp5=0.0, temp6 = 0.0;
3974+ double temp7 = 0.0, temp8=0.0, temp9 = 0.0;
3975+ for(int k=0; k<OrbitalType_end; k++){
3976+ for(int l=0; l<OrbitalType_end; l++){
3977+
3978+ temp1 += tmpRotMat2ndDerivs [i][k][dimA1][dimA2]
3979+ *tmpRotMat [j][l]
3980+ *tmpDiaOverlapAOsInDiaFrame[k][l];
3981+ temp2 += tmpRotMat [i][k]
3982+ *tmpRotMat2ndDerivs [j][l][dimA1][dimA2]
3983+ *tmpDiaOverlapAOsInDiaFrame[k][l];
3984+ temp3 += tmpRotMat [i][k]
3985+ *tmpRotMat [j][l]
3986+ *tmpDiaOverlapAOs2ndDerivs [k][l][dimA1][dimA2];
3987+ temp4 += tmpRotMat1stDerivs [i][k][dimA1]
3988+ *tmpRotMat1stDerivs [j][l][dimA2]
3989+ *tmpDiaOverlapAOsInDiaFrame[k][l];
3990+ temp5 += tmpRotMat1stDerivs [i][k][dimA1]
3991+ *tmpRotMat [j][l]
3992+ *tmpDiaOverlapAOs1stDerivs [k][l][dimA2];
3993+ temp6 += tmpRotMat1stDerivs [i][k][dimA2]
3994+ *tmpRotMat1stDerivs [j][l][dimA1]
3995+ *tmpDiaOverlapAOsInDiaFrame[k][l];
3996+ temp7 += tmpRotMat [i][k]
3997+ *tmpRotMat1stDerivs [j][l][dimA1]
3998+ *tmpDiaOverlapAOs1stDerivs [k][l][dimA2];
3999+ temp8 += tmpRotMat1stDerivs [i][k][dimA2]
4000+ *tmpRotMat [j][l]
4001+ *tmpDiaOverlapAOs1stDerivs [k][l][dimA1];
4002+ temp9 += tmpRotMat [i][k]
4003+ *tmpRotMat1stDerivs [j][l][dimA2]
4004+ *tmpDiaOverlapAOs1stDerivs [k][l][dimA1];
40164005 }
4017-
4018- diatomicOverlapAOs2ndDerivs[i][j][dimA1][dimA2] = temp1+temp2+temp3
4019- +temp4+temp5+temp6
4020- +temp7+temp8+temp9;
40214006 }
4007+
4008+ diatomicOverlapAOs2ndDerivs[i][j][dimA1][dimA2] = temp1+temp2+temp3
4009+ +temp4+temp5+temp6
4010+ +temp7+temp8+temp9;
40224011 }
40234012 }
40244013 }
4025-
40264014 }
4027- catch(MolDSException ex){
4028- this->FreeDiatomicOverlapAOs2ndDeriTemps(&diaOverlapAOsInDiaFrame,
4029- &diaOverlapAOs1stDerivInDiaFrame,
4030- &diaOverlapAOs2ndDerivInDiaFrame,
4031- &rotMat,
4032- &rotMat1stDerivatives,
4033- &rotMat2ndDerivatives,
4034- &tempDiaOverlapAOs1stDerivs,
4035- &tempDiaOverlapAOs2ndDerivs);
4036- throw ex;
4037- }
4038- // free
4039- this->FreeDiatomicOverlapAOs2ndDeriTemps(&diaOverlapAOsInDiaFrame,
4040- &diaOverlapAOs1stDerivInDiaFrame,
4041- &diaOverlapAOs2ndDerivInDiaFrame,
4042- &rotMat,
4043- &rotMat1stDerivatives,
4044- &rotMat2ndDerivatives,
4045- &tempDiaOverlapAOs1stDerivs,
4046- &tempDiaOverlapAOs2ndDerivs);
4015+
40474016 /*
40484017 for(int i=0; i<OrbitalType_end; i++){
40494018 for(int j=0; j<OrbitalType_end; j++){
@@ -4058,9 +4027,25 @@ void Cndo2::CalcDiatomicOverlapAOs2ndDerivatives(double**** diatomicOverlapAOs2n
40584027 }
40594028
40604029 void Cndo2::CalcDiatomicOverlapAOs2ndDerivatives(double**** diatomicOverlapAOs2ndDerivs,
4030+ double** tmpDiaOverlapAOsInDiaFrame,
4031+ double** tmpDiaOverlapAOs1stDerivInDiaFrame,
4032+ double** tmpDiaOverlapAOs2ndDerivInDiaFrame,
4033+ double*** tmpDiaOverlapAOs1stDerivs,
4034+ double**** tmpDiaOverlapAOs2ndDerivs,
4035+ double** tmpRotMat,
4036+ double*** tmpRotMat1stDerivs,
4037+ double**** tmpRotMat2ndDerivs,
40614038 int indexAtomA,
40624039 int indexAtomB) const{
40634040 this->CalcDiatomicOverlapAOs2ndDerivatives(diatomicOverlapAOs2ndDerivs,
4041+ tmpDiaOverlapAOsInDiaFrame,
4042+ tmpDiaOverlapAOs1stDerivInDiaFrame,
4043+ tmpDiaOverlapAOs2ndDerivInDiaFrame,
4044+ tmpDiaOverlapAOs1stDerivs,
4045+ tmpDiaOverlapAOs2ndDerivs,
4046+ tmpRotMat,
4047+ tmpRotMat1stDerivs,
4048+ tmpRotMat2ndDerivs,
40644049 *this->molecule->GetAtom(indexAtomA),
40654050 *this->molecule->GetAtom(indexAtomB));
40664051 }
@@ -4084,42 +4069,6 @@ double Cndo2::Get2ndDerivativeElementFromDistanceDerivatives(double firstDistanc
40844069 return value;
40854070 }
40864071
4087-void Cndo2::MallocDiatomicOverlapAOs2ndDeriTemps(double*** diaOverlapAOsInDiaFrame,
4088- double*** diaOverlapAOs1stDerivInDiaFrame,
4089- double*** diaOverlapAOs2ndDerivInDiaFrame,
4090- double*** rotMat,
4091- double**** rotMat1stDerivs,
4092- double***** rotMat2ndDerivs,
4093- double**** tempDiaOverlapAOs1stDerivs,
4094- double***** tempDiaOverlapAOs2ndDerivs) const{
4095- MallocerFreer::GetInstance()->Malloc<double>(diaOverlapAOsInDiaFrame, OrbitalType_end, OrbitalType_end);
4096- MallocerFreer::GetInstance()->Malloc<double>(diaOverlapAOs1stDerivInDiaFrame, OrbitalType_end, OrbitalType_end);
4097- MallocerFreer::GetInstance()->Malloc<double>(diaOverlapAOs2ndDerivInDiaFrame, OrbitalType_end, OrbitalType_end);
4098- MallocerFreer::GetInstance()->Malloc<double>(rotMat, OrbitalType_end, OrbitalType_end);
4099- MallocerFreer::GetInstance()->Malloc<double>(rotMat1stDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end);
4100- MallocerFreer::GetInstance()->Malloc<double>(rotMat2ndDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end, CartesianType_end);
4101- MallocerFreer::GetInstance()->Malloc<double>(tempDiaOverlapAOs1stDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end);
4102- MallocerFreer::GetInstance()->Malloc<double>(tempDiaOverlapAOs2ndDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end, CartesianType_end);
4103-}
4104-
4105-void Cndo2::FreeDiatomicOverlapAOs2ndDeriTemps(double*** diaOverlapAOsInDiaFrame,
4106- double*** diaOverlapAOs1stDerivInDiaFrame,
4107- double*** diaOverlapAOs2ndDerivInDiaFrame,
4108- double*** rotMat,
4109- double**** rotMat1stDerivs,
4110- double***** rotMat2ndDerivs,
4111- double**** tempDiaOverlapAOs1stDerivs,
4112- double***** tempDiaOverlapAOs2ndDerivs) const{
4113- MallocerFreer::GetInstance()->Free<double>(diaOverlapAOsInDiaFrame, OrbitalType_end, OrbitalType_end);
4114- MallocerFreer::GetInstance()->Free<double>(diaOverlapAOs1stDerivInDiaFrame, OrbitalType_end, OrbitalType_end);
4115- MallocerFreer::GetInstance()->Free<double>(diaOverlapAOs2ndDerivInDiaFrame, OrbitalType_end, OrbitalType_end);
4116- MallocerFreer::GetInstance()->Free<double>(rotMat, OrbitalType_end, OrbitalType_end);
4117- MallocerFreer::GetInstance()->Free<double>(rotMat1stDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end);
4118- MallocerFreer::GetInstance()->Free<double>(rotMat2ndDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end, CartesianType_end);
4119- MallocerFreer::GetInstance()->Free<double>(tempDiaOverlapAOs1stDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end);
4120- MallocerFreer::GetInstance()->Free<double>(tempDiaOverlapAOs2ndDerivs, OrbitalType_end, OrbitalType_end, CartesianType_end, CartesianType_end);
4121-}
4122-
41234072 // calculate OverlapAOs matrix. E.g. S_{\mu\nu} in (3.74) in J. A. Pople book by GTO expansion.
41244073 // See Eqs. (28) - (32) in [DY_1977]
41254074 void Cndo2::CalcOverlapAOsByGTOExpansion(double** overlapAOs,
--- a/src/cndo/Cndo2.h
+++ b/src/cndo/Cndo2.h
@@ -215,9 +215,25 @@ protected:
215215 int indexAtomA,
216216 int indexAtomB) const;
217217 void CalcDiatomicOverlapAOs2ndDerivatives(double**** overlapAOs2ndDeri,
218+ double** tmpDiaOverlapAOsInDiaFrame,
219+ double** tmpDiaOverlapAOs1stDerivInDiaFrame,
220+ double** tmpDiaOverlapAOs2ndDerivInDiaFrame,
221+ double*** tmpDiaOverlapAOs1stDerivs,
222+ double**** tmpDiaOverlapAOs2ndDerivs,
223+ double** tmpRotMat,
224+ double*** tmpRotMat1stDerivs,
225+ double**** tmpRotMat2ndDerivs,
218226 const MolDS_base_atoms::Atom& atomA,
219227 const MolDS_base_atoms::Atom& atomB) const;
220228 void CalcDiatomicOverlapAOs2ndDerivatives(double**** overlapAOs2ndDeri,
229+ double** tmpDiaOverlapAOsInDiaFrame,
230+ double** tmpDiaOverlapAOs1stDerivInDiaFrame,
231+ double** tmpDiaOverlapAOs2ndDerivInDiaFrame,
232+ double*** tmpDiaOverlapAOs1stDerivs,
233+ double**** tmpDiaOverlapAOs2ndDerivs,
234+ double** tmpRotMat,
235+ double*** tmpRotMat1stDerivs,
236+ double**** tmpRotMat2ndDerivs,
221237 int indexAtomA,
222238 int indexAtomB) const;
223239 double Get2ndDerivativeElementFromDistanceDerivatives(double firstDistanceDeri,
@@ -496,25 +512,8 @@ private:
496512 double**** diisStoredErrorVect,
497513 double*** diisErrorProducts,
498514 double** diisErrorCoefficients);
499- void MallocDiatomicOverlapAOs2ndDeriTemps(double*** diaOverlapAOsInDiaFrame,
500- double*** diaOverlapAOs1stDerivInDiaFrame,
501- double*** diaOverlapAOs2ndDerivInDiaFrame,
502- double*** rotMat,
503- double**** rotMat1stDerivs,
504- double***** rotMat2ndDerivs,
505- double**** tempDiaOverlapAOs1stDerivs,
506- double***** tempDiaOverlapAOs2ndDerivs) const;
507- void FreeDiatomicOverlapAOs2ndDeriTemps(double*** diaOverlapAOsInDiaFrame,
508- double*** diaOverlapAOs1stDerivInDiaFrame,
509- double*** diaOverlapAOs2ndDerivInDiaFrame,
510- double*** rotMat,
511- double**** rotMat1stDerivs,
512- double***** rotMat2ndDerivs,
513- double**** tempDiaOverlapAOs1stDerivs,
514- double***** tempDiaOverlapAOs2ndDerivs) const;
515515 };
516516
517-
518517 }
519518 #endif
520519
--- a/src/mndo/Mndo.cpp
+++ b/src/mndo/Mndo.cpp
@@ -981,6 +981,9 @@ void Mndo::MallocTempMatricesEachThreadCalcHessianSCF(double***** diatomicOve
981981 double****** tmpDiatomicTwoElecTwoCore1stDerivs,
982982 double*** tmpDiaOverlapAOsInDiaFrame,
983983 double*** tmpDiaOverlapAOs1stDerivInDiaFrame,
984+ double*** tmpDiaOverlapAOs2ndDerivInDiaFrame,
985+ double**** tmpDiaOverlapAOs1stDerivs,
986+ double***** tmpDiaOverlapAOs2ndDerivs,
984987 double*** tmpRotatedDiatomicOverlap,
985988 double*** tmpMatrix) const{
986989 MallocerFreer::GetInstance()->Malloc<double>(diatomicOverlapAOs1stDerivs,
@@ -1041,6 +1044,18 @@ void Mndo::MallocTempMatricesEachThreadCalcHessianSCF(double***** diatomicOve
10411044 MallocerFreer::GetInstance()->Malloc<double>(tmpDiaOverlapAOs1stDerivInDiaFrame,
10421045 OrbitalType_end,
10431046 OrbitalType_end);
1047+ MallocerFreer::GetInstance()->Malloc<double>(tmpDiaOverlapAOs2ndDerivInDiaFrame,
1048+ OrbitalType_end,
1049+ OrbitalType_end);
1050+ MallocerFreer::GetInstance()->Malloc<double>(tmpDiaOverlapAOs1stDerivs,
1051+ OrbitalType_end,
1052+ OrbitalType_end,
1053+ CartesianType_end);
1054+ MallocerFreer::GetInstance()->Malloc<double>(tmpDiaOverlapAOs2ndDerivs,
1055+ OrbitalType_end,
1056+ OrbitalType_end,
1057+ CartesianType_end,
1058+ CartesianType_end);
10441059 MallocerFreer::GetInstance()->Malloc<double>(tmpRotatedDiatomicOverlap,
10451060 OrbitalType_end, OrbitalType_end);
10461061 MallocerFreer::GetInstance()->Malloc<double>(tmpMatrix,
@@ -1060,6 +1075,9 @@ void Mndo::FreeTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverl
10601075 double****** tmpDiatomicTwoElecTwoCore1stDerivs,
10611076 double*** tmpDiaOverlapAOsInDiaFrame,
10621077 double*** tmpDiaOverlapAOs1stDerivInDiaFrame,
1078+ double*** tmpDiaOverlapAOs2ndDerivInDiaFrame,
1079+ double**** tmpDiaOverlapAOs1stDerivs,
1080+ double***** tmpDiaOverlapAOs2ndDerivs,
10631081 double*** tmpRotatedDiatomicOverlap,
10641082 double*** tmpMatrix) const{
10651083 MallocerFreer::GetInstance()->Free<double>(diatomicOverlapAOs1stDerivs,
@@ -1120,6 +1138,18 @@ void Mndo::FreeTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverl
11201138 MallocerFreer::GetInstance()->Free<double>(tmpDiaOverlapAOs1stDerivInDiaFrame,
11211139 OrbitalType_end,
11221140 OrbitalType_end);
1141+ MallocerFreer::GetInstance()->Free<double>(tmpDiaOverlapAOs2ndDerivInDiaFrame,
1142+ OrbitalType_end,
1143+ OrbitalType_end);
1144+ MallocerFreer::GetInstance()->Free<double>(tmpDiaOverlapAOs1stDerivs,
1145+ OrbitalType_end,
1146+ OrbitalType_end,
1147+ CartesianType_end);
1148+ MallocerFreer::GetInstance()->Free<double>(tmpDiaOverlapAOs2ndDerivs,
1149+ OrbitalType_end,
1150+ OrbitalType_end,
1151+ CartesianType_end,
1152+ CartesianType_end);
11231153 MallocerFreer::GetInstance()->Free<double>(tmpRotatedDiatomicOverlap,
11241154 OrbitalType_end,
11251155 OrbitalType_end);
@@ -1678,9 +1708,13 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
16781708 double***** tmpDiatomicTwoElecTwoCore1stDerivs = NULL;
16791709 double** tmpDiaOverlapAOsInDiaFrame = NULL; // diatomic overlapAOs in diatomic frame
16801710 double** tmpDiaOverlapAOs1stDerivInDiaFrame = NULL; // first derivative of the diaOverlapAOs. This derivative is related to the distance between two atoms.
1711+ double** tmpDiaOverlapAOs2ndDerivInDiaFrame = NULL; // second derivative of the diaOverlapAOs. This derivative is related to the distance between two atoms.
1712+ double*** tmpDiaOverlapAOs1stDerivs = NULL; // first derivatives of the diaOverlapAOs. This derivatives are related to the all Cartesian coordinates.
1713+ double**** tmpDiaOverlapAOs2ndDerivs = NULL; //sedond derivatives of the diaOverlapAOs. This derivatives are related to the all Cartesian coordinates.
16811714 double** tmpRotMat1stDeriv = NULL;
16821715 double** tmpRotatedDiatomicOverlap = NULL;
16831716 double** tmpMatrix = NULL;
1717+
16841718 try{
16851719 this->MallocTempMatricesEachThreadCalcHessianSCF(&diatomicOverlapAOs1stDerivs,
16861720 &diatomicOverlapAOs2ndDerivs,
@@ -1694,6 +1728,9 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
16941728 &tmpDiatomicTwoElecTwoCore1stDerivs,
16951729 &tmpDiaOverlapAOsInDiaFrame,
16961730 &tmpDiaOverlapAOs1stDerivInDiaFrame,
1731+ &tmpDiaOverlapAOs2ndDerivInDiaFrame,
1732+ &tmpDiaOverlapAOs1stDerivs,
1733+ &tmpDiaOverlapAOs2ndDerivs,
16971734 &tmpRotatedDiatomicOverlap,
16981735 &tmpMatrix);
16991736 #pragma omp for schedule(auto)
@@ -1717,6 +1754,14 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
17171754 indexAtomA,
17181755 indexAtomB);
17191756 this->CalcDiatomicOverlapAOs2ndDerivatives(diatomicOverlapAOs2ndDerivs[indexAtomB],
1757+ tmpDiaOverlapAOsInDiaFrame,
1758+ tmpDiaOverlapAOs1stDerivInDiaFrame,
1759+ tmpDiaOverlapAOs2ndDerivInDiaFrame,
1760+ tmpDiaOverlapAOs1stDerivs,
1761+ tmpDiaOverlapAOs2ndDerivs,
1762+ tmpRotMat,
1763+ tmpRotMat1stDerivs,
1764+ tmpRotMat2ndDerivs,
17201765 indexAtomA,
17211766 indexAtomB);
17221767 this->CalcDiatomicTwoElecTwoCore1stDerivatives(diatomicTwoElecTwoCore1stDerivs[indexAtomB],
@@ -1798,6 +1843,9 @@ void Mndo::CalcHessianSCF(double** hessianSCF, bool isMassWeighted) const{
17981843 &tmpDiatomicTwoElecTwoCore1stDerivs,
17991844 &tmpDiaOverlapAOsInDiaFrame,
18001845 &tmpDiaOverlapAOs1stDerivInDiaFrame,
1846+ &tmpDiaOverlapAOs2ndDerivInDiaFrame,
1847+ &tmpDiaOverlapAOs1stDerivs,
1848+ &tmpDiaOverlapAOs2ndDerivs,
18011849 &tmpRotatedDiatomicOverlap,
18021850 &tmpMatrix);
18031851 }// end of omp-region
--- a/src/mndo/Mndo.h
+++ b/src/mndo/Mndo.h
@@ -159,6 +159,9 @@ private:
159159 double****** tmpDiatomicTwoElecTwoCore1stDerivs,
160160 double*** tmpDiaOverlapAOsInDiaFrame,
161161 double*** tmpDiaOverlapAOs1stDerivInDiaFrame,
162+ double*** tmpDiaOverlapAOs2ndDerivInDiaFrame,
163+ double**** tmpDiaOverlapAOs1stDerivs,
164+ double***** tmpDiaOverlapAOs2ndDerivs,
162165 double*** tmpRotatedDiatomicOverlap,
163166 double*** tmpMatrix) const;
164167 void FreeTempMatricesEachThreadCalcHessianSCF(double***** diatomicOverlapAOs1stDerivs,
@@ -173,6 +176,9 @@ private:
173176 double****** tmpDiatomicTwoElecTwoCore1stDerivs,
174177 double*** tmpDiaOverlapAOsInDiaFrame,
175178 double*** tmpDiaOverlapAOs1stDerivInDiaFrame,
179+ double*** tmpDiaOverlapAOs2ndDerivInDiaFrame,
180+ double**** tmpDiaOverlapAOs1stDerivs,
181+ double***** tmpDiaOverlapAOs2ndDerivs,
176182 double*** tmpRotatedDiatomicOverlap,
177183 double*** tmpMatrix) const;
178184 double GetAuxiliaryHessianElement1(int mu,