CHARMM Development Project Forums CHARMM Interfaces QM/MM Discussion and Questions SCCDFTB not converged in qm/mm minimization step

Dear All,

Currently I am working on a system where two ligands and one amino acid residue are treated quantum mechanically. The following error was obtained when the minimization step ran.

***** LEVEL -5 WARNING FROM *****
***** SCCDFTB not converged.
******************************************
BOMLEV ( -1) IS REACHED - TERMINATING. WRNLEV IS 5

The respective .inp and .out file are attached. What does it mean from this error and how to get rid off this?.

Thanks.

Wild uninformed guess, is the net charge correct? Or an unusual/incorrect protonation state that the TB method is struggling with? I haven't used this method with CHARMM, hence no clue about expected output, but it's a bit strange it doesn't seem to get through even one SCF iteration.

There was a mistake in the input file (definition of the R1and H1) .So I made it correct and ran again.But the same error was obtained.I have no idea where it goes wrong since I am new to qm/mm.I saw some similar post where it says the SHAKE is causing this problem and as a solution, it was said to turn off the shake off (reduce the timestep). But in this script, I guess it is not possible.If anyone has used sccdftb method with QM/MM minimization can you please please guide me in order to correct this error?

Input file
*SCCDFTB
*
bomblev -2 !For parameters
! Topology file
! -------------
stream toppar.str
stream addlink.str

open read unit 18 card name "step1_sys_tcts_with_link_atoms.psf" !changing psf to cmp1
read psf unit 18 card

set RADIUS 40

set sys tcts_sccdftb ! system from labbook ??


open unit 15 read card name "step1_sys_tcts_with_link_atoms.crd"
read coor card unit 15
close unit 15

define QM sele (segid R1 .or. segid H1) show end !
define MM sele .not. QM end
!endif

! delete atoms sele resname TIP3 end

cons fix sele none end
!mini abnr nstep 100


open write form unit 13 name check.pdb
write coor pdb unit 13
close unit 13

scalar WMAIN set 1.0 sele (QM) .and. type P* SHOW end ! main coordinate weights
scalar WMAIN set 2.0 sele (QM) .and. type O* SHOW end
scalar WMAIN set 3.0 sele (QM) .and. type N* SHOW end
scalar WMAIN set 4.0 sele (QM) .and. type C* SHOW end
scalar WMAIN set 5.0 sele (QM) .and. type H* SHOW end
scalar WMAIN set 5.0 sele (QM) .and. type QQH1 SHOW end
!endif
! Do not forget about the link atom

update atom fshift noextend cdie vdw vfswitch eps 1.0 - ! constant dielectric and classical force shift method
cutnb 14 ctofnb 12 ctonnb 10.0 -
WMIN 1.5 inbfrq -1


!WE WILL NOW GO THROUGH A SERIES OF MINIMIZATIONS

define backbone sele segid RIC .and. (type C .or. type CA .or. type O .or. type N) end
define include sele .byres. (point 0.0 0.0 0.0 cut @RADIUS) end
!Keep outside fix ALL the time
cons fix sele .not. include end ! Everything outside is fixed

shake bonh tol 1.0e-10 param sele MM end


!NOW WE WILL WORK WITH INCLUDE

SCCDFTB REMOve SELE QM end - ! classical energies within QM atoms are removed
TEMP 0.00 SCFT 0.000001 EXGR - ! electronic temp (fermi) , convergence criteria
chrg -2.00 HBON D3RD mull ! MIXE 1 !slash NOC RC 50 BA 2 LT 5 !without dispersion

!#######################################################

!#################################################
energy

open read unit 4 formatted name wat45.pot
sbound read unit 4
close unit 4

sbound set xref 0.0 yref 0.0 zref 0.0 assign 1 -
sele (resname TIP3 .and. type OH2) end


!SETUP REACTION COORDINATES
!#########################################################
!reaction distance


!###################################################
cons fix sele none end
cons fix sele .not. include end
mini abnr nstep 100 tolegr 0.05

cons fix sele .not. include end

open write form unit 13 name "@SYS.crd"
write coord card unit 13
* Minimized QMMM
*

write coord pdb name "@SYS.pdb" sele all end

open write form unit 15 name "@SYS.psf"
write psf card unit 15
stop

the part of the output file where the error was obtained

CHARMM>

CHARMM> !#######################################################
CHARMM>

CHARMM> !#################################################
CHARMM> energy

NONBOND OPTION FLAGS:
ELEC VDW ATOMs CDIElec FSHIft VATOm VFSWIt
BYGRoup NOEXtnd NOEWald
CUTNB = 14.000 CTEXNB =999.000 CTONNB = 10.000 CTOFNB = 12.000
CGONNB = 0.000 CGOFNB = 10.000
WMIN = 1.500 WRNMXD = 0.500 E14FAC = 1.000 EPS = 1.000
NBXMOD = 5
There are 13220460 atom pairs and 70627 atom exclusions.
There are 0 group pairs and 8566 group exclusions.
Generating nonbond list with Exclusion mode = 5
== PRIMARY == SPACE FOR 17724330 ATOM PAIRS AND 0 GROUP PAIRS

General atom nonbond list generation found:
13192711 ATOM PAIRS WERE FOUND FOR ATOM LIST
783507 GROUP PAIRS REQUIRED ATOM SEARCHES

ewevge: ier = 155

***** LEVEL -5 WARNING FROM *****
***** SCCDFTB not converged.
******************************************
BOMLEV ( -1) IS REACHED - TERMINATING. WRNLEV IS 5

Again I'm backing a structural reason. Is FAD the Flavin Adenine Dinucleotide coenzyme? If so, it's a pretty floppy molecule so if it's nowhere near a minimum, a TB method probably won't get you there. Have you done any MM minimisation prior to this? TB methods struggle with species that have unusual protonation states or are unusually charged so it may be that it just don't deal with your phosphates very well.

If you're sure of your geometry, maybe apply the classic hardware maintenance logic; unplug everything except the bare minimum. So try it with only a couple of QM atoms and see if you can get it to converge then add more to isolate the problem.