Kimball Model of EthylCyanoAcetate

We will write a simple Mathematica version 10 code to study Ethylcyanoacetate.Our object is to compute its energy components given a 3D structure from the PDB library and to plot the different electronic constituents.
All computations are transparent and operations explained in detail. This is an HTML copy of the original Mathematica Notebook EtCyanoAc_n1.nb, which can be downloaded and interactively used.
(ES 27 March 2013/26 Sept. 2014)

Input and Definitions

The stored coordinates are read from the Wolfram Mathematica ChemicalData repository (coordinates in pm = picometer) into the t1 array. We use atomic units here, the universally applied system of theoretical chemistry and (micro)physics, see NIST. Length data are in Bohr: 1 a0 = 0.52917721 Å = 52.917721 pm (whose reciprocal is 0.01889726 bohr/pm, used below), electric charges in ± electron charges, and energies in Hartrees: 1 Eh = 2 Rydberg = 627.5095 kcal/mol = 2’625.50 kJ/mol. Other units are indicated explicitely.

EtCyanoAc_nexp_1.gif

EtCyanoAc_nexp_2.gif

EtCyanoAc_nexp_3.gif

EtCyanoAc_nexp_4.gif

EtCyanoAc_nexp_5.gif

EtCyanoAc_nexp_6.gif

EtCyanoAc_nexp_7.gif

Sort array of atoms, show charges, number of “heavy” atoms, number of hydrogens

EtCyanoAc_nexp_8.gif

EtCyanoAc_nexp_9.gif

EtCyanoAc_nexp_10.gif

EtCyanoAc_nexp_11.gif

EtCyanoAc_nexp_12.gif

Compute Kimball radii from distance matrix

EtCyanoAc_nexp_13.gif

EtCyanoAc_nexp_14.gif

Select sigma bonded pairs by a distance criterion

EtCyanoAc_nexp_15.gif

EtCyanoAc_nexp_16.gif

Build vector of core radii and proton excentricities from CH4, NH3, H2O gauge molecules (cnofhydb.pas)

EtCyanoAc_nexp_17.gif

EtCyanoAc_nexp_18.gif

From lower triangle of symmetrical matrix of distances subtract vector of core radii

EtCyanoAc_nexp_19.gif

EtCyanoAc_nexp_20.gif

EtCyanoAc_nexp_21.gif

Transpose to upper triangular matrix, show, then subtract vector of core radii from the other side, which leaves diameter of bonding clouds between “heavy” atoms and radius of H-clouds:

EtCyanoAc_nexp_22.gif

EtCyanoAc_nexp_23.gif

EtCyanoAc_nexp_24.gif

EtCyanoAc_nexp_25.gif

Determine total number of σ-bonds and radii of bonding clouds between “heavy” atoms

EtCyanoAc_nexp_26.gif

EtCyanoAc_nexp_27.gif

EtCyanoAc_nexp_28.gif

Form some arrays, to be used later

EtCyanoAc_nexp_29.gif

Compute kinetic energy of sigma bonding and core clouds. Kimballs ansatz 1.125/R2 per electron

EtCyanoAc_nexp_31.gif

EtCyanoAc_nexp_32.gif

EtCyanoAc_nexp_33.gif

EtCyanoAc_nexp_34.gif

Total kinetic energy except for π-clouds and lone pairs

EtCyanoAc_nexp_35.gif

EtCyanoAc_nexp_36.gif

EtCyanoAc_nexp_37.gif

EtCyanoAc_nexp_38.gif

EtCyanoAc_nexp_39.gif

Determine connectivity matrix:

EtCyanoAc_nexp_40.gif

EtCyanoAc_nexp_41.gif

EtCyanoAc_nexp_42.gif

Localize double bonds and positions of π-clouds:

Transform the triangle of every target atom with two of its neighbors into the xy-plane and attach π-clouds above and below the plane to the target. Then backtransform the triangle plus π-clouds into the molecular coordinate array. This is done by the two subroutines PItrans.m, and PICNtrans.m, that can be downloaded from the URL’s given below for reading them in.

EtCyanoAc_nexp_43.gif

Localize lone pairs, compute orientation:

Transform the triangle of every target atom with two of its neighbors into the xy-plane and attach lone pair(s). Then backtransform the triangle plus lone pair(s) into molecular coordinate array. Some such transformations shown below for perpendicular and in-plane oxygen lone pairs. One side of the triangle is parallel to x-axis, the target atom at origin. Subroutines: NCNtrans.m, ODtrans.m, and OBWtrans.m.

EtCyanoAc_nexp_44.gif

EtCyanoAc_nexp_45.gif

EtCyanoAc_nexp_46.gif

EtCyanoAc_nexp_47.gif

Bonding clouds: Atom pair, radius of cloud; statistics of sizes

EtCyanoAc_nexp_48.gif

EtCyanoAc_nexp_49.gif

EtCyanoAc_nexp_50.gif

EtCyanoAc_nexp_51.gif

Plot of EthylCyanoacetate and its partial constituents

EtCyanoAc_nexp_52.gif

Graphics:EthylCyanoacetate

Graphics:Core skeleton

Graphics:bonded σ skeleton

Graphics:π - clouds on skeleton

Graphics:Lone Pairs

Preparation of interaction matrices (w-w, w-n, check size arrays):

EtCyanoAc_nexp_58.gif

EtCyanoAc_nexp_59.gif

EtCyanoAc_nexp_60.gif

EtCyanoAc_nexp_61.gif

EtCyanoAc_nexp_62.gif

Compute energy components (check charges):

EtCyanoAc_nexp_63.gif

EtCyanoAc_nexp_64.gif

EtCyanoAc_nexp_65.gif

Politzerratio

EtCyanoAc_nexp_66.gif

EtCyanoAc_nexp_67.gif

Results (energies in [Eh] Hartree)

EtCyanoAc_nexp_68.gif

EtCyanoAc_nexp_69.gif

EtCyanoAc_nexp_70.gif

EtCyanoAc_nexp_71.gif

EtCyanoAc_nexp_72.gif

EtCyanoAc_nexp_73.gif

EtCyanoAc_nexp_74.gif

EtCyanoAc_nexp_75.gif

EtCyanoAc_nexp_76.gif

EtCyanoAc_nexp_77.gif

EtCyanoAc_nexp_78.gif

EtCyanoAc_nexp_79.gif

EtCyanoAc_nexp_80.gif

EtCyanoAc_nexp_81.gif

EtCyanoAc_nexp_82.gif

EtCyanoAc_nexp_83.gif

EtCyanoAc_nexp_84.gif

Created with the Wolfram Language