[molpro-user] DFT calculations produced different energy with a different guesses
Kirk Peterson
kipeters at wsu.edu
Wed Nov 3 21:26:15 GMT 2010
David,
you could try as well to do the HF calculation with the multi program to see if this has better convergence characteristics. I often do this for difficult cases. Just specify the same occ and closed as your final occupations below. I agree that it's strange your orbitals are not ordered by energy. Are you sure it converged correctly?
-Kirk
On Nov 3, 2010, at 1:45 PM, David Danovich wrote:
> Dear Kirk,
>
> In general, you are right and it is always possible to converge to different electronic state, especially for the systems with transition metals. For DFT methods I always do TDDFT calculations in order to be sure that the state I get is a lower one (for such systems with 3 transition metal atoms I have no ideas how to do CASSCF with 33 valence electrons and I do not know a priory which d electrons are important). But for this particular case I think I have the same state with the same occupancy as you can see below and also with almost the same eigenvectors. But as I said before difference in total energy is 3.8 kcal/mol. What is also strange that for irreducible representation A1 Aufbau principle is not satisfy. You may see that orbital 11.1 has much lower energy than orbitals 6.1, 7.1, 8.1, 9.1, 10.1 for both calculations. My feeling is that something goes wrong with DFT calculations in Molpro. I will try to carry out this type of calculations with other programs.
>
> Thank you David
>
> Case 1)
>
> Final alpha occupancy: 12 5 9 4
> Final beta occupancy: 10 5 8 4
>
> Orbital energies:
>
> 1.1 2.1 3.1 4.1 5.1 6.1 7.1 8.1 9.1 10.1
> -4.671414 -4.118564 -3.065611 -2.328939 -2.323338 -0.438461 -0.379638 -0.343993 -0.313299 -0.288287
>
> 11.1 12.1 13.1 14.1
> -0.509678 -0.218071 -0.055982 -0.024797
>
> 1.2 2.2 3.2 4.2 5.2 6.2 7.2
> -3.028087 -2.319030 -0.478500 -0.324411 -0.306541 -0.066386 0.004508
>
> 1.3 2.3 3.3 4.3 5.3 6.3 7.3 8.3 9.3 10.3
> -4.118395 -3.028437 -2.326183 -2.321552 -0.481322 -0.301593 -0.288912 -0.282174 -0.196122 -0.041215
>
> 11.3
> 0.035134
>
> 1.4 2.4 3.4 4.4 5.4 6.4
> -2.318707 -0.440376 -0.298955 -0.284142 0.018513 0.139495
>
> Case 3)
>
> Final alpha occupancy: 12 5 9 4
> Final beta occupancy: 10 5 8 4
> 1.1 2.1 3.1 4.1 5.1 6.1 7.1 8.1 9.1 10.1
> -4.673062 -4.117760 -3.014878 -2.327371 -2.323105 -0.441349 -0.360882 -0.345206 -0.319350 -0.288610
>
> 11.1 12.1 13.1 14.1
> -0.514110 -0.216705 -0.051875 -0.024126
>
> 1.2 2.2 3.2 4.2 5.2 6.2 7.2
> -3.061378 -2.318440 -0.460149 -0.323223 -0.306348 -0.070494 0.001935
>
> 1.3 2.3 3.3 4.3 5.3 6.3 7.3 8.3 9.3 10.3
> -4.117591 -3.050251 -2.324891 -2.321034 -0.451631 -0.300020 -0.288446 -0.281090 -0.195381 -0.042282
>
> 11.3
> 0.035544
>
> 1.4 2.4 3.4 4.4 5.4 6.4
> -2.318118 -0.492930 -0.298691 -0.283211 0.018322 0.139423
>
>
>
>
> On Wed, Nov 3, 2010 at 9:40 PM, Kirk Peterson <kipeters at wsu.edu> wrote:
> David,
>
> it is certainly not always the case that the resulting HF or DFT solution is independent of the initial guess. In a perfect world that would be nice, but often a poorly constructed initial guess (or just bad luck) can result in convergence to a local minimum, i.e., an upper electronic state. This happens quite often for transition metal complexes because of the plethora of possible d-orbital occupations. The only way to get a feeling for what is going on is to print the orbitals and inspect the higher lying occupied ones. Is the bonding as you might expect? Another option is to carry out a CASSCF calculation as I suggested before in order to determine what the character of the higher lying occupied orbitals should really be. One can always then use the rotate command to swap in the correct orbitals in your HF or DFT calculation.
>
> -Kirk
>
> On Nov 3, 2010, at 9:03 AM, David Danovich wrote:
>
>> Dear Kirk,
>>
>> You are right, the calculations converge to different structures. But the question is why. In all calculations as you can see in the input I sent in my original post I used the same starting geometry, just guess wave functions for B3LYP 4B2 state were different. I got different energy already in the first step (at the point where geometry was the same for all three cases). It means that even for single pint calculations solution is depending on starting guess. In general, I can use any guess wave function and should get the same converged result. But as you see I got different results. Occupancy and eigenvalues are more or less similar for the cases 1) and 3) but in any case the difference in energy is around 3.8 kcal/mol. How one can be sure that he really get the lower solution (for particular geometry)? Do you have any suggestion.
>>
>> Thank you David
>>
>>
>> On Wed, Nov 3, 2010 at 5:01 PM, Kirk Peterson <kipeters at wsu.edu> wrote:
>> David,
>>
>> it looks to me that you're getting stuck into a few local minima. Do these different calculations converge to different structures? How different are the energies at the starting geometries? I would strongly recommend printing the orbitals and taking a careful look. You could also try doing a CASSCF calculation (print the orbitals and CI vector) to see how things should really be occupied.
>>
>> regards,
>>
>> Kirk
>>
>>
>> On Nov 3, 2010, at 6:19 AM, David Danovich wrote:
>>
>>> Hello,
>>>
>>> I am calculating CuAu2 molecule in high spin state (4B2 state) using B3LYP method. Below you can find three different possibilities I have used in the calculations. In the first one 1) I have calculated first HF wave function and used it as a guess for DFT calculation. Energy I got was -468.86467917 au. In the second calculation 2) I directly have calculated B3LYP without HF wave function. The energy I got was -468.83995914 au. In the third calculation 3) I first have done B3LYP calculation for state (4A2) and then used it as a guess for the calculations of 4B2. The energy I got was
>>> -468.85824088 au. As you can see there is quite large difference in energy for the same state. I was trying to use different grids but it does not solved the problem. In my opinion result should not depend on the guess so drastically. What can be a solution for this problem?
>>>
>>> Thank you in advance David
>>>
>>> 1)
>>>
>>> ***, Peterson PP tz basis set
>>> memory,250,m
>>> r = 2.72684681 ang;
>>> a = 60.65125012 degree;
>>> geometry={Cu1; !z-matrix geometry input
>>> Au2,Cu1,r;
>>> Au3,Cu1,r,Au2,a;
>>> }
>>> basis=cc-pVTZ-PP
>>>
>>> {hf,maxit=500;
>>> wf,57,3,3}
>>>
>>> {rks,b3lyp3,maxit=500;
>>> wf,57,3,3}
>>>
>>> optg
>>>
>>> Total energy -468.86467917 au
>>> ____________________________________________________
>>>
>>> 2)
>>> ***, Peterson PP tz basis set
>>> memory,250,m
>>> r = 2.72684681 ang;
>>> a = 60.65125012 degree;
>>> geometry={Cu1; !z-matrix geometry input
>>> Au2,Cu1,r;
>>> Au3,Cu1,r,Au2,a;
>>> }
>>> basis=cc-pVTZ-PP
>>>
>>> {ks,b3lyp3,maxit=500;
>>> shift,-0.1,0.0;
>>> wf,57,3,3}
>>>
>>> optg
>>>
>>> Total energy -468.83995914 au
>>>
>>> ______________________________________________________
>>>
>>> 3)
>>> ***, Peterson PP tz basis set
>>> memory,250,m
>>> r = 2.72684681 ang;
>>> a = 60.65125012 degree;
>>> geometry={Cu1; !z-matrix geometry input
>>> Au2,Cu1,r;
>>> Au3,Cu1,r,Au2,a;
>>> }
>>>
>>> ! cc-pVTZ-PP
>>>
>>> basis=cc-pVTZ-PP
>>>
>>> {ks,b3lyp3,maxit=500;
>>> wf,57,4,3}
>>>
>>> {ks,b3lyp3,maxit=500;
>>> shift,-1.0,0.0;
>>> wf,57,3,3}
>>>
>>> optg
>>>
>>> Total energy -468.85824088 au
>>>
>>> _______________________________________________
>>> Molpro-user mailing list
>>> Molpro-user at molpro.net
>>> http://www.molpro.net/mailman/listinfo/molpro-user
>>
>>
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.molpro.net/pipermail/molpro-user/attachments/20101103/8a061073/attachment.html>
More information about the Molpro-user
mailing list