Showing posts with label internal coordinates. Show all posts
Showing posts with label internal coordinates. Show all posts

05 September 2013

512. Briefly: zmatrices in nwchem -- methanol

And another update:
I can now confirm that using your own z matrix still does not constrain the geometry during a PES scan, which was the original impetus for this post: http://verahill.blogspot.com.au/2013/09/511-when-nwchem-pes-scans-fail-to.html

Another update:
the gaussian run failed after 14 geometry steps during the first PES point.
NTrRot= -1 NTRed= 628 NAtoms= 34 NSkip= 532 IsLin=F Error in internal coordinate system. Error termination via Lnk1e in /opt/gaussian/g09/l103.exe at Thu Sep 5 18:17:12 2013. Job cpu time: 0 days 22 hours 25 minutes 27.6 seconds. File lengths (MBytes): RWF= 192 Int= 0 D2E= 0 Chk= 28 Scr= 1
Not being an expert, to me it seems that there's something fundamentally difficult with the system I'm working on. In an ideal world I'd give the actual details, but quite apart from the risk of being scooped, doing so would also make it easier to identify me (not that it's impossible at this point).

[Suffice to say that the system holds a large polyoxoanion and a small p-block anion, both of which are symmetrical and negatively charged. The goal of the PES scan is to bring the ions closer to see whether they 'react'. Which is also a troublesome use of computational resources -- computational chemistry is good at answering well-defined questions using carefully designed computational experiments -- but not generally very good at answering ill-defined questions about synthesis (i.e. you can't generally 'mix two things together and see what happens' and expect a useful result. Anyway, regardless of that, that's exactly what I want to do.]

Update:
nwchem still gives errors about autoz in spite of using noautoz. But I also get messages about the user generated z matrix, so we'll see whether my input is respected or not.

Also, for one of the calcs I'm getting
There are insufficient internal variables: expected 95 got 96

which is really, really, really annoying since there doesn't seem to be a real fix for it -- I've tried everything suggested in http://www.nwchem-sw.org/index.php/Special:AWCforum/st/id286. I can get the same calc to run in gaussian though (gaussian has its own issues), but it'd be nice if stuff just...worked...

Original post:
Normally you don't have to fiddle with zmatrices in nwchem -- instead you'd typically supply cartesian coordinates, and nwchem would do autoz to autogenerate internal (z matrix) coordinates.

Sometimes that fails, and nwchem defaults to using cartesian coordinates. In most cases, this isn't a cause for any real concern -- the computation will continue although I think cartesian coordinates are supposed to be slightly slower.

However, if you're doing a PES scan you'll notice that it's not proceeding as intended -- the constraints are completely ignored: 511. When nwchem PES scans fail to constrain -- autoz failure

The easiest remedy is to supply the internal coordinates directly, but there honestly aren't too many examples online showing how that's done, and I kept on getting annoying failure messages along the lines of
NWChem Input Module ------------------- zmat ---- THE 3-D PIECE OF -Z- DATA FOR ATOM = 2 IS NEITHER FLOATING POINT NOR ALPHANUMERIC OR COULD NOT BE MATCHED WITH A VARIABLE. STOP IAT= 2 ZMAT= 2 1 0 0 0 0.00000 0.00000 0.00000 ------------------------------------------------------------------------ JOB STOPPED PROGRAM STOP IN - ZDAT - ------------------------------------------------------------------------ ------------------------------------------------------------------------ CALLS IT QUIT FROM HND_HNDERR 0 ------------------------------------------------------------------------ This error has not yet been assigned to a category

This particular error came about because the zmatrix module is case sensitive, and my Variables couldn't be interpreted (it should be variables). Anyway, you'll understand more after this post, and it isn't important anyway.


Calculation using a z matrix (internal coordinates) in nwchem, with a little bit of help from openbabel:

Assuming that you set up a calculation in e.g. ECCE for a geometry optimisation of methanol you'll end up with the following input file:
scratch_dir /home/andy/scratch Title "methanol" Start methanol echo charge 0 geometry autosym units angstrom C 0.0351714 0.00548884 0.0351714 H -0.617781 -0.634073 0.667983 H 0.667983 -0.634073 -0.617781 H -0.605139 0.646470 -0.605139 O 0.839603 0.818768 0.839603 H 1.38912 0.201564 1.38912 end ecce_print ecce.out basis "ao basis" spherical print H library "6-31+G*" O library "6-31+G*" C library "6-31+G*" END dft mult 1 direct XC b3lyp grid fine mulliken end driver default end task dft optimize
Take the coordinates, and paste them into a file, e.g. methanol.xyz:
6 methanol C 0.0351714 0.00548884 0.0351714 H -0.617781 -0.634073 0.667983 H 0.667983 -0.634073 -0.617781 H -0.605139 0.646470 -0.605139 O 0.839603 0.818768 0.839603 H 1.38912 0.201564 1.38912

Next, use openbabel:
babel -ixyz methanol.xyz -ogzmat 
#Put Keywords Here, check Charge and Multiplicity. methanol 0 1 C H 1 r2 H 1 r3 2 a3 H 1 r4 2 a4 3 d4 O 1 r5 2 a5 3 d5 H 5 r6 1 a6 2 d6 Variables: r2= 1.1117 r3= 1.1117 a3= 109.74 r4= 1.1094 a4= 108.78 d4= 118.90 r5= 1.3984 a5= 110.18 d5= 238.51 r6= 0.9924 a6= 105.98 d6= 60.61 1 molecule converted 18 audit log messages
The format isn't quite right (everything in red needs to go, and the V in blue should be lower case), but we can sort that out:

 babel -ixyz ~/methanol.xyz -ogzmat |sed 's/\=//g;s/V/v/g;s/\://g' |tail -n+6 > methanol.zmat
C H 1 r2 H 1 r3 2 a3 H 1 r4 2 a4 3 d4 O 1 r5 2 a5 3 d5 H 5 r6 1 a6 2 d6 variables r2 1.1117 r3 1.1117 a3 109.74 r4 1.1094 a4 108.78 d4 118.90 r5 1.3984 a5 110.18 d5 238.51 r6 0.9924 a6 105.98 d6 60.61

Let's update out nwchem input file with the internal coordinates:
scratch_dir /home/andy/scratch Title "methanol" Start methanol echo charge 0 geometry noautoz zmatrix C H 1 r2 H 1 r3 2 a3 H 1 r4 2 a4 3 d4 O 1 r5 2 a5 3 d5 H 5 r6 1 a6 2 d6 variables r2 1.1117 r3 1.1117 a3 109.74 r4 1.1094 a4 108.78 d4 118.90 r5 1.3984 a5 110.18 d5 238.51 r6 0.9924 a6 105.98 d6 60.61 end end ecce_print ecce.out basis "ao basis" spherical print H library "6-31+G*" O library "6-31+G*" C library "6-31+G*" END dft mult 1 direct XC b3lyp grid fine mulliken end driver default end task dft optimize

And run. Done!