As always, Uwe comes and saves the day
as usual, there are more things possible with Turbomole than one might expect
Which is hard to find, given the scarce information on the documentation...
1. yes, you can fix cartesian coordinates and optimize in redundant coordinates at the same time. This is an undocumented feature, and I am not sure if it is stable in your case and how efficient it will be, but if you fix the cartesian coords, define internal redundant coordinates (with ired in define), and then run the job with jobex -statpt [...], i.e. using the new statpt module instead of relax, it might work. At least with Turbomole 5.10
Ok, we have turbomole 5.10 here so I tried this. I do not see any information anywhere on the fixed coordinates (not if I go into the 'i' menu in define or if I check the coordinates file). If I do it by hand, through 'idef' I see some fixed coordinates in a separate keyword ($intdef) in the coordinates file. So I may be doing something wrong. If this works as you say, then it's the thing I want. Didn't seem to work for me though...
2. z-matrix for Turbomole? molden, babel, gaussview (via xyz format), or even the tmole script (which is part of Turbomole, see documentation chapter 1.9) can do that. But please remember that a z-matrix is a very intuitive method to define coordinates, but mathematically not really proper (one has to add dummy atoms to be able to describe a set of cartesian coordinates in some cases). Hence, using internal (redundant) coordinates is much more efficient than using z-matrices.
Yes, I agree and I usually use internal redundant coordinates. But, if the description by z-matrices saves me a huge amount of time and allows me to do it more automatically, then in a few cases it may be worth it.
3. what is the reason for your SCF convergence problems? My experience with very difficult cases is that it might take a while to find the correct occupation and DIIS settings, etc. But once you have converged orbitals, a geometry optimization if often less demanding than the first energy calculation. Especially if only the hydrogens are allowed to move. Did you try $fermi ?
I haven't tried $fermi yet. I've been working with other software packages (which I am more familiar with) trying to get this calculation going but to no avail. So, I'm on my first few tries with Turbomole... the SCF problem comes from my experience with the other programs, not especially from the results of Turbomole usage. I'll give a look at $fermi though.
4. depending on which method you are using, the coordinates of experimental crystal structures do not have to be a good start structure for the optimization. The minimum geometry of HF, DFT, MP2, etc. can differ significantly from X-ray structures.
I know and ultimately the idea is to optimize the whole thing. But since we are going to compare the results with some solid stare NMR calculations and NMR is highly dependent on the geometry I want to use something as similar to what they have as possible. I'll also optimize the full geometry and see if it gets anything the experimentalists were not able to explain.
5. did you try to optimize the hole structure without constrains? RI-DFT with B-P86 functional, $marij and def-SV(P) basis set should not take too long for 94 atoms - but it depends if you have transition metal atoms, etc., so you might first try to get converged orbitals (see point 3).
I do not have transition metals (I have some Na
+ though) and I did try that. I increased the SCF steps to 200 and the calculation stop after 200 SCF steps without converging...
As I said before, I went through the idef and fixed stuff by hand. I wanted to see the behavior of the program if I fixed stuff...
Since my idea is to fix anything that has (only) the atoms 1 through 60, it is just a matter of getting to the file and seeing which lines have only those numbers, so it shouldn't be hard to script.
All goes well if I look only and bond lengths. However, for angles and dihedrals the thing screws up a little bit... Here go the two sets of bonds, angles and dihedrals I fixed.
$intdef
# definitions of internal coordinates
1 f 1.0000000000000 stre 18 25 val= 2.80071
2 f 1.0000000000000 bend 14 16 13 val= 106.00700
3 f 1.0000000000000 tors 10 14 13 12 val= -13.82794
4 f 1.0000000000000 stre 6 5 val= 2.52608
5 f 1.0000000000000 bend 3 1 2 val= 125.73701
6 f 1.0000000000000 tors 45 44 41 40 val= 1.93098
And here is the same for the part automatically made by define (I'll change the numbering on the first column to the same numbers as above, just to make it more clear):
1 k 1.0000000000000 stre 18 25 val= 2.80071
2 k -0.6161601823403 bend 14 16 13 val= 20.30057
3 k -0.8090169943750 tors 10 14 13 12 val= 25.24627
3 k -0.5877852522925 tors 10 14 13 12 val= -18.09457
4 k 1.0000000000000 stre 6 5 val= 2.52608
5 k -0.5000000000000 bend 3 1 2 val= 3.16414
5 k 0.8660254037844 bend 3 1 2 val= 4.79577
5 k -1.0000000000000 bend 3 1 2 val= -4.46234
6 k -0.5000000000000 tors 45 44 41 40 val= -8.27260
6 k 0.8660254037844 tors 45 44 41 40 val= -12.92675
6 k 1.0000000000000 tors 45 44 41 40 val= 8.54507
So, as you see, the 'val' is totally different! The manual says the internal coordinates do not have a direct meaning and cannot be changed by the user (page 46, first paragraph) but I do not want to change them, I want to fix them! Is there a way to get from the bottom information to the top one? Of course if I get your first suggestion to work, this is totally unnecessary...
Comments from anyone?
Thank you very much for your help, Uwe.
Best regards,
Daniel