1 / 72

Utilizing AICD and GIMIC programs to study magnetically induced current density

Utilizing AICD and GIMIC programs to study magnetically induced current density. Sobereva 2012-Jul-26. “Left-hand rule“ to determine the direction of the currents induced by a magnetic field in a conductor. dia=Diatropic para=Paratropic.

emanion
Download Presentation

Utilizing AICD and GIMIC programs to study magnetically induced current density

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript

  1. Utilizing AICD and GIMIC programs to study magnetically induced current density Sobereva 2012-Jul-26

  2. “Left-hand rule“ to determine the direction of the currents induced by a magnetic field in a conductor

  3. dia=Diatropic para=Paratropic

  4. Aromatic molecules: The diatropic component on the outside of the molecules dominates over the paratropic one inside the ring yielding a net diatropic current • Antiaromatic molecules: The paratropic current inside the ring dominates • Nonaromatic molecules: Also sustain diatropic and paratropic ring currents, but the two components cancel, yielding a minuscule net current

  5. Maps of the current density provide qualitative information about the current pathways and the calculated current strengths complete the quantitative picture of the current delocalization. Porphin 1 bohr above Magnetic field is perpendicular to the molecular plane

  6. Studying induced current density in quantitative manner Current strength: Numerical integration of the current density passing specific bonds Phys. Chem. Chem. Phys., 2011, 13, 20500–20518

  7. Broader applications of current strength E: Energy of hydrogen bond JHB: Strength of the current passing hydrogen bond Current strength can be used as a measure of the electron delocalization thus indicating how well the investigated molecule is able to transport charges

  8. Visualizing induced current density in different way Anisotropy of the Induced Current Density (AICD) • Advantages of AICD relative to current density map • A scalar field • Linearly independent of the total electron density • Independent of the relative orientation of the molecule and the magnetic field • Not only for aromatic systems but also for any kind of conjugation (through bond, through space, ...) in any kind of system (ground state, excited state, transition state, ...)

  9. Limitation of current density analysis: Some delocalization regions are invisible Cancelled out

  10. Magnetically induced current density T: Induced current density tensor (matrix) J and T dependent on real space coordinate

  11. Alpha current density Beta current density Total current density Spin current density

  12. Problem in calculation of magnetical properties: Gauge origin Momentum operator: In magnetic field, generalized form: (q=-1 for electron) A is vector potential due to magnetic field Gauge transformation: (Curl of any gradient of scalar potential is zero) equivalently, RO is arbitrarily selected Gauge origin

  13. For complete basis-set, magnetical properties (NMR shielding, induced current density, etc.) are uniquely determined by external magnetic field, while for finite basis-set, gauge origin problem should be explicitly taken into account Solutions: • GIAO (Gauge-Including Atomic Orbital) • IGLO (Individual Gauge for Localized Orbitals) • CSGT (Continuous Set of Gauge Transformations) • IGAIM (Individual Gauge for Atoms In Molecules) GIAO is the default method for NMR calculation in most popular quantum chemistry programs, the result converges faster with respect to basis-set than CSGT

  14. The way to generate induced current density • Anisotropy of the Induced Current Density (AICD) program Using T directly, which is an intermediate quantity in CSGT calculation • Gauge including magnetically induced current (GIMIC) program , component of NMR shielding tensor for nucleus K m: nuclear magnetic moment Commonly, σis calculated by derivative technique, in which coupled perturbed Hartree–Fock (CPHF) method is used to generate derivative of density matrix

  15. : Spin ,,: Cartesian components ,: Basis function index h: One-electron Hamiltonian : GIAO (which eliminated gauge origin dependency)

  16. AICD program Features: 1. Calculating AICD in specific spatial region 2. Generating POVRAY input file, by which isosurface of AICD and induced current density vectors on the isosurface can be visualized 3. Orbital contributions can be decomposed (sigma-pi separation can be realized)

  17. Written in C++ • Very difficult to be hacked, a sample: • Minimax Grenzen (for_each (Stromdichte.begin (), Stromdichte.end (), Minimax ())) ; • Matrix<float> Mittelpunkt (Null,1,3) ; • Mittelpunkt = ((*Stromdichte.begin()).GetOrt ()+(*(++Stromdichte.begin())).GetOrt ())/ float (2.0) ; • float Radius = 0 ; • { • for (TStromdichten::iterator Index = Stromdichte.begin(); Index != Stromdichte.end(); ++Index) • Radius = max(Radius,(Mittelpunkt - (*Index).GetOrt ()).Laenge ()) ; • } • Free, open source, can be obtained from the author: rherges@oc.uni-kiel.de • Require T exported by NMR=CSGT calculation in Gaussian98/03. However T is unable to be outputted by public version, source code of G98/03 is needed to be patched and then recompiled. • Latest version: 1.5.7.1

  18. Compile AICD program • Uncompress AICD-1.5.7.1.tar.bz2 • Enter the folder, run make all • ln -s ./AICD /usr/local/bin/AICD GNU compiler is enough

  19. Patching Gaussian 98/03 • Compile Gaussian03 source code (one can consult http://hi.baidu.com/sobereva/blog/item/abe4d744fec12b85b3b7dc52.html) • Copy l1002.F.g03.diff from AICD-1.5.7.1/gaussian_patch to g03 folder • Backup original l1002.F and l1002.exe • Enter g03 folder and run patch < l1002.F.g03.diff • Search “Call FFRead(IEOF)” in l1002.F, delete the line “goto 12” behind it. Search “If(V(IOccMO+I-1).eq.Two) then”, add below two lines after the line “Write(IOut,1010) I” • else • Write(IOut,1020) I • 6. Enter csh environment by typing csh, enter g03 folder, run source bsd/g03.login, then run mg l1002.exe to compile new l1002.exe

  20. Usage Print full help information: AICD –longhelp Typical command: AICD -p 80000 -pov -b 0 0 1 -m 4 benzene_AICD.out Most useful options -p 80000: Number of grid points -pov: Output input file of POVRAY render -b 0 0 1: Magnetic field vector, will be normalized automatically -m 2: Single view point -m 4: Multiple view points -l 0.05: Isovalue of AICD -s: Smooth isosurface of AICD (time consuming) -rot 90 0 0: Rotate view angle

  21. Example: Benzene #P b3lyp/6-31G* NMR=csgt b3lyp/6-31+g(d) opted 0 1 C 0.00000000 1.39864200 0.00000000 C 1.21125900 0.69932100 0.00000000 C 1.21125900 -0.69932100 0.00000000 C 0.00000000 -1.39864200 0.00000000 C -1.21125900 -0.69932100 0.00000000 C -1.21125900 0.69932100 0.00000000 H 0.00000000 2.48606800 0.00000000 H 2.15299800 1.24303400 0.00000000 H 2.15299800 -1.24303400 0.00000000 H 0.00000000 -2.48606800 0.00000000 H -2.15299800 -1.24303400 0.00000000 H -2.15299800 1.24303400 0.00000000 17 17 20 21 B3LYP/6-31G* is often enough to yield qualitatively satisfactory result Only pi orbitals (MO17, MO20 and MO21) will be taken into accounted. By default all orbitals will be considered

  22. Run: AICD -p 80000 -pov -b 0 0 1 -m 4 benzene_AICD.out • Output file: • benzene_AICD.icd.gz: Compressed version of benzene_AICD.icd, which recorded T information • benzene_AICD.icd80000.gz: Binary version of AICD cube file • benzene_AICD.icd80000.cub: cube file of AICD, can be visualized by VMD, Multiwfn, Molekel, etc. • benzene_AICD.pdb: Molecule geometry file • POVRAY input file: • benzene_AICD_80000_0.050_0_0_1_Aniso_4.4.Molekuel.inc • benzene_AICD_80000_0.050_0_0_1_Aniso_4.4.Isoober.inc • benzene_AICD_80000_0.050_0_0_1_Aniso_4.4.Rotate.inc • benzene_AICD_80000_0.050_0_0_1_Aniso_4.4.Molekuel.pov • benzene_AICD_80000_0.050_0_0_1_Aniso_4.4.RenderMich.pov

  23. Visualizing result 1. Download Windows version of POVRAY and install it 2. Copy the five POVRAY input files to a new folder 3. Double-click the file icon of benzene_AICD_80000_0.050_0_0_1_Aniso_4.4.RenderMich.pov to load it into POVRAY 4. Adjust rendering setting in top-left of POVRAY window, then click Run button POVRAY is the most popular render, can be freely download at http://www.povray.org/download/

  24. Taking all orbitals into account AICD -p 80000 -pov -b 0 0 1 -m 2 -s -rot 45 0 0 benzene_AICD_all.out

  25. For more details about AICD, see Chem. Rev., 105, 3758 http://hi.baidu.com/sobereva/blog/item/88b6fa544c64e43e3a293581.html

  26. GIMIC program • Functions: • Outputting induced current density, its modulus, divergence in specific 2D plane or 3D spatial region • Outputting electron density in specific 2D plane or 3D spatial region • Integration of the current flow through defined cut-planes in molecules

  27. Written in Fortran90, easy to be hacked • Free, open source, can be downloaded at https://gitorious.org/gimic • Principally can be interfaced to any quantum chemistry programs that with GIAO implmented. Currently only interfaced to CFOUR/ACESII and Turbomole • Unable to decompose orbital contribution like AICD

  28. Digression: A brief introduction of CFOUR Main page: http://www.cfour.de Manual: http://slater.chemie.uni-mainz.de/cfour/index.php?n=Main.Manual • CFOUR=Coupled-Cluster techniques for Computational Chemistry • Renamed from ACESII-MAB, which is a variant of ACESII • Free, open source, but fax or regular mail is mandatory to apply a license • High efficient coupled-cluster calculation, analytic derivatives are available for many sophisticated theoretical methods • Interface to GIMIC, MRCC, DIRAC, NEWTON-X • Lacking of a well written manual

  29. Available Hamiltonians R/U/ROHF TCSCF MP2/3/4 CC2/3 CCD CCSD CCSD+T CCSD(T) CCSDT CCSDT-n (n=1-4) Brueckner-CCD (B-CCD) and B-CCD(T) Orbital-optimized CC (OO-CC) CID/CISD QCISD and QCISD(T) DFT is not supported!

  30. Analytic gradient HF-SCF TCSCF MP2 MP3 MP4 CC2 CCD QCISD CCSD QCISD(T) CCSD(T) CCSDT-n (n=1-4) CC3 CCSDT Analytic second derivative HF-SCF (RHF, UHF, ROHF) MP2 (RHF and UHF) MP3 (RHF and UHF) MP4 (RHF and UHF) CCD (RHF and UHF) CCSD (RHF and UHF) CCSD(T) (RHF and UHF) CCSDT-n (n=1-4, RHF) CCSDT (RHF) Task: Single point, geometry optimization, locating TS, excitation state, etc. Properties: IR, Raman, anharmonic, NMR, g-tensor, etc.

  31. Compilation of serial version of CFOUR Decompress cfour_v1.tar.gz and enter the folder ./configure FC=ifort --enable-gimic (or FC=gfortran) make Then executable files can be found in cfour_v1/bin folder, add the path to $PATH environment variable, e.g. export PATH=$PATH:/sob/cfour_v1/bin

  32. Compilation of parallel version of CFOUR with MKL Install intel fortran compiler with MKL to default location (version 10.x, 11.x,12.x are OK) Install mpich2 to default location (version 1.4.1p1 is OK) export MKLPATH=/opt/intel/mkl/lib/intel64 ./configure FC=ifort MPIFC=mpif90 --enable-gimic --with-blas="$MKLPATH/libmkl_solver_ilp64.a -Wl,--start-group $MKLPATH/libmkl_intel_ilp64.a $MKLPATH/libmkl_intel_thread.a $MKLPATH/libmkl_core.a -Wl,--end-group -openmp -lpthread" --enable-mpi=mpich --with-mpirun="mpirun -np \$CFOUR_NUM_CORES" --with-exenodes="mpirun -np \$CFOUR_NUM_CORES“ If ifort 12.x is used, change “-O3” in make.config to “-O2” make

  33. Example of input file Water CC-LR/DZP at experimental equilibrium geometry O H 1 R H 1 R 2 A R=0.958 A=104.5 *CFOUR(CALC=CCSD,BASIS=DZP,EXCITE=EOMCC) %excite* 1 1 1 5 0 6 0 1.0 [Blank line]

  34. Input file should be named as ZMAT GENBAS and ECPDATA (if pseudopotential is employed) should be present in the same folder as ZMAT, the two files can be found in cfour_v1/basis For parallel implementation, CC_PROGRAM=ECC and ABCDTYPE=AOBASIS must be specified in ZMAT, and set below environment variables export CFOUR_NUM_CORES=4 export MKL_NUM_THREADS=1 In the folder containing ZMAT, run: xcfour |tee xxx.out

  35. Full key words list: http://slater.chemie.uni-mainz.de/cfour/index.php?n=Main.FileStructure MEMORY=xxx将可用内存设为xxx，默认单位为INTEGERWORDS。对于32/64bit平台分别乘以4/8就是KB。用MEM_UNIT可以将单位改为kB, MB, GB, and TB。 ECP=ON：使用赝势。 CHARGE：体系的电荷 MULTIPLICITY：体系多重度 COORDINATES：控制体系的坐标描述。默认的INTERNAL是内坐标，CARTESIAN可以用笛卡尔坐标，但是不能做几何优化。XYZ2INT是提供z坐标连接关系但使用笛卡尔坐标描述位置。 SCF_CONV=N：密度矩阵最大变化小于10^-N就停了。默认为7。 SCF_DAMPING：设500有益于解决SCF不收敛 SCF_MAXCYC：SCF最大迭代次数，默认150 SCF_EXTRAPOLATION：是否用DIIS，默认为ON SPHERICAL：默认的ON是用球谐型高斯函数。OFF用笛卡尔型 SUBGROUP：默认使用最高阶的阿贝尔点群对称性。C1就相当于SYMMETRY=OFF

  36. 只要在相应变量上写上星号，就代表对它们进行优化，其它的关键词不需要。只要在相应变量上写上星号，就代表对它们进行优化，其它的关键词不需要。 O H 1 R* H 1 R* 2 A* GEO_METHOD=TS：找过渡态。 GEO_CONV=N：设定优化收敛限为N Hartree/bohr。 GEO_MAXCYC：最大优化步数，默认为50。 GEO_MAXSTEP：设几何优化的最大步长为millbohr。默认300。 EVAL_HESS=N：每隔N个优化步算一次精确Hessian，默认为从不，也就是用准牛顿法。 EXCITATION=n：计算第n激发态。 SCF_CONV：控制SCF收敛。 ABCDTYPE=AOBASIS：建议对最高至CCSD(T)的各种计算都加上，可以加快速度。默认是=0 (STANDARD) XFIELD,YFIELD,ZFIELD：XYZ方向加电场 PRINT=1：可以比默认0少输出些信息 PROPS=FIRST_ORDER：计算一阶属性（多极矩，相对论校正，电场梯度，自旋密度，Mulliken电荷）；SECOND_ORDER计算静态可极化率；DYNAMICAL；计算含频可极化率；NMR：计算NMR；HYPERPOL：计算静态超极化率；DYN_HYP计算含频率超极化率 RAMAN_INT=1：算拉曼 RELATIVISTIC：设定相对论校正 VIBRATION=ANALYTIC：用解析二阶导数计算谐振频率

  37. CC_PROGRAM：控制做CC的程序，默认是VCC，对于CCSD, CCSD+T, CCSD(T), closed-shell CCSDT-n, CC3, and CCSDT，建议用ECC来加快速度。 DBOC=1：做diagonal Born-Oppenheimer对能量的校正。 FROZEN_CORE：CFOUR默认在post-HF中不冻芯。如果设为ON，则内核轨道在post-HF过程中都被冻结。如果想设定具体哪些被冻，则用DROPMO来设定。 DROPMO：如果输入比如1>10-55-58>64，就代表1,2,3,4,5,6,7,8,9,10,55,58,59,60,61,62,63轨道都被冻结。 FROZEN_VIRT：默认为OFF，如果设为ON，高于指定能量的虚轨道将不被考虑。 HFSTABILITY：ON代表做SCF波函数稳定性测试 REFERENCE：RHF、UHF、ROHF、TCSCF(two-configureational SCF) RESTART_CC=1：重启CC计算。 CC_CONV=N：CC收敛标准为系数最大改变值小于10^-N。默认为7。

  38. Modifying GIMIC code In cubeplot.f90, search do i=1,npts(1) do j=1,npts(2) do k=1,npts(3) write(fd,'(f12.6)',advance='no') pdata(i,j,k) if (mod(l,6) == 5) write(fd,*) l=l+1 end do end do end do Replace them by do i=1,npts(1) do j=1,npts(2) write(fd,"(6(1PE13.5))",advance="no") pdata(i,j,1:npts(3)) write(fd,*) end do end do Change 0 to 1 in line 35 Add write(fd,"(i5,4f12.6)") 1,1D0,0D0,0D0,0D0 behind line 38

  39. Add below code behind line 610 in jfield.f if (mod(l,6) /= 5) then l=0 write(fd1,*) write(fd2,*) end if Search “write(fd, '(i5,3f12.6)') 0, origin”, change 0 to 1 Search “write(fd, '(i5,3f12.6)') npts(3), 0.d0, 0.d0, step(3)”, add write(fd,"(i5,4f12.6)") 1,1D0,0D0,0D0,0D0 behind it

  40. Compile and install GIMIC ./configure FC=ifort make gimic.x will be produced in GIMIC folder Open GIMIC/tools/MOL2mol.sh, modify ACES2 varible to actual path, e.g. ACES2=/sob/cfour_v1. Change two $ACES in the file to $ACES2. Open GIMIC/tools/xgimic2.sh, modify ACES2 varible to actual path make install Copy GIMIC/tools/xgimic2.sh to /usr/local/bin

  41. Necessary files for running GIMIC • mol: Molecular geometry and basis-set • XDENS: Density matrix and its derivative with respect to magnetic field • gimic.inp: GIMIC input file

  42. Typical CFOUR input file for GIMIC BENZENE X C 1 RCC C 1 RCC 2 A60 C 1 RCC 3 A60 2 D180 C 1 RCC 4 A60 3 D180 C 1 RCC 5 A60 4 D180 C 1 RCC 6 A60 5 D180 H 1 RXH 2 A60 7 D180 H 1 RXH 3 A60 2 D180 H 1 RXH 4 A60 3 D180 H 1 RXH 5 A60 4 D180 H 1 RXH 6 A60 5 D180 H 1 RXH 7 A60 6 D180 RCC=1.39864 RXH=2.48607 A60=60.0 D180=180.0 *CFOUR(CALC=HF,BASIS=6-31G*,PROPS=NMR,SYMMETRY=ON) [Blank line]

  43. 1. In the folder containing the ZMAT file, run xgimic2.sh --scf |tee cfour.out (Use --cc and –mbpt instead if couple-cluster and MP2/3/4 is used, respectively) This script invokes CFOUR to generate XDENS 2. Run MOL2mol.sh to produce mol file 3. Prepare gimic.inp based on template, which can be found in GIMIC/examples 4. Finally, run gimic

  44. Keywords in GIMIC input file Location of GIMIC manual: GIMIC/Documentation/gimic.pdf All length units are Bohr dryrun=off # don't actually calculate (good for tuning grids, etc.) mpirun=off # run in parallel mode title="" basis="mol" # Name of MOL file with coordinates and basis sets density="XDENS" # File with AO density and perturbed densities spherical=off # don't touch, unless you REALLY know what you are doing debug=1 # debug print level diamag=on # turn on/off diamagnetic contributions paramag=on # turn on/off paramagnetic contributions GIAO=on # turn on/off GIAOs. Don't change unless you know why. openshell=false screening=on # use screening to speed up screen_thrs=1.d-8 # Screening threshold show_up_axis=true # mark "up" axis in .xyz files calc=[cdens, integral, divj] ……

  45. Induced current density cdens { magnet_axis=-z #[-] i,j,k || x,y,z -> align magnet along axis # magnet=[0.0, -1.0, 0.0] # magnet vector scale_vectors=1.0 jtensor="JTENSOR" # file name for the current tensors jvector="JVECTOR" # file name for the current vectors grid(base) { # grid type can be: base/std, bond or file (see below) type=even # even spaced or gauss quadrature distribution on interval origin=[-8.0, -8.0, 0] # origin of grid ivec=[1.0, 0.0, 0.0] # basis vector i jvec=[ 0.0, 1.0, 0.0] # basis vector j ( k = i x j ) lengths=[16.0, 16.0, 16.0] # lenthts of (i,j,k) spacing=[0.2, 0.2, 0.5] # spacing of points on grid (i,j,k) # grid_points=[50,50,0] # number of gridpoints on grid (i,j,k) rotation=[0.0,0.0,0.0] # Rotation of (i,j,k) -> (i',j',k') in degrees # Euler angles, x-y-z convention } plot { # file names for plots vector="JVEC" modulus="JMOD" cube_mod="jmod" # gopenmol="jmod.plt" #projection="JPRJ" } }

  46. Integration current flow integral { # magnet_axis=T magnet=[0.0, 0.0, 1.0] modulus=off # calculate the |J| integral spin=total # total, alpha, beta grid(bond) { # define grid orthogonal to a bond type=gauss # gauss distribution of grid points # origin=[0.0, 0.0, 0.0] # fix grid orientation # bond=[1,2] # atom indeces for bond # fixpoint=5 bond=[1,2] # atom indeces for bond fixpoint=5 # coord1=[0.0, 0.0, 2.145166] # "atom" coordinates # coord2=[0.0, 0.0, -2.145166] # fixcoord=[0.0, 0.0, 0.0] distance=1.323 # place grid 'distance' between atoms ## Either grid_points or spacing can be specified. The number of grid ## points will be rounded upwards to nearest multiple of gauss_order gauss_order=9 # order for gauss quadrature grid_points=[30,30,0] # number of gridpoints on grid (i,j,k) # spacing=[0.1, 0.1, 0.1] # spacing of points on grid (i,j,k) up=6.0 down=6.0 in=2.0 out=6.0 # height=[-5.0, 5.0] # width=[-5.0, 5.0] rotation=[0.0,0.0,0.0] # Rotation of (i,j,k) -> (i',j',k') # radius=1.0 # round grid: cut off at radius } }

More Related