gdirect gthresh,energy=1d-8,orbital=1d-7,grid=1d-8 symmetry,nosym angstrom geometry={ 1,O1,, -1.551007, -0.114520, 0.000000 2,H2,, -1.934259, 0.762503, 0.000000 3,H3,, -0.599677, 0.040712, 0.000000 4,O4,, 1.350625, 0.111469, 0.000000 5,H5,, 1.680398, -0.373741, -0.758561 6,H6,, 1.680398, -0.373741, 0.758561} basis={ set,orbital; default,avdz set,jkfit; default,avdz/jkfit set,mp2fit; default,avdz/mp2fit} oldnorm=1 ca=2110.2; cb=2120.2 nkerna=1; nkernb=1 {df-hf,basis=jkfit; save,2100.2} edm_hf=energy dummy,4,5,6 {df-hf,basis=jkfit; save,$ca; start,atdens} ema_hf=energy monomerA,orb=$ca,core=$nkerna dummy,1,2,3 {df-hf,basis=jkfit; save,$cb; start,atdens} emb_hf=energy monomerB,orb=$cb,core=$nkernb {grid,gridthr=1d-6} {dfsapt,e1p=1,e1x=1,e2i=1,e2xi=1,auxbas='jkfit',fexx=1.0,fcoul=1.0,mode=1,cdirect=1,thrcg=1d-6} {grid,gridthr=1d-8} eint_hf=(edm_hf-ema_hf-emb_hf)*1000e0 mH delta_hf=eint_hf-e12tot dummy,4,5,6 lxgamma=1.15 {df-ks,lxbecke,pbex,pbec,basis=jkfit; dftfac,0.25,0.75,1.0; asymp,0d0; save,$ca; start,atdens} monomerA,orb=$ca,core=$nkerna dummy,1,2,3 lxgamma=1.15 {df-ks,lxbecke,pbex,pbec,basis=jkfit; dftfac,0.25,0.75,1.0; asymp,0d0; save,$cb; start,atdens} monomerB,orb=$cb,core=$nkernb {grid,gridthr=1d-6} {dfsapt,all=1,auxbas='mp2fit',fxc=1.0,fcoul=1.0,mode=1,cpks=2,fxcfit=1,fxcsing=1,t678=3} eint_sapt=e12tot+delta_hf