How to build 2Fe-2S cluster model for simulation?

User discussions
,
1

GROMACS version: 2022.5
Related software: ACPYPE, Ambertools22, Gaussian16

The protein I am interested in is PDB:6B6G. There is a 2Fe-2S cluster locating at the center of the protein and the two ion atoms interact with two deprotonated CYS repectively.

I already get some parameters from literature and DFT calculation related to the 2Fe2S - 4CYS. Now I would like to embed these bond parameter into ff99sb-ildn force field. How should I build this model for the gromacs simulation? Which of the topology files should be modified?

Thanks for helping!

2

I write several .itp files. First is the topol.top:
; Include forcefield parameters
#include “amber99sb.ff/forcefield.itp”
; Include custom atom type
#include “dummy.itp”
#include “atomtype_FES.atp”
; Include chain topologies
#include “topol_chain_C_modified.itp”
#include “topol_chain_D_modified.itp”
; Include FES cluster topology
#include “FeS.itp”
#ifdef POSRES_FES
#include “posre_FeS.itp”
#endif
#include “CYS-Fes.itp”

For dummy.itp, I just replace the H atom of Cys residue by dummy atom. It is shown here:
[ atomtypes ]
; name at.num mass charge ptype sigma epsilon
DUM 0 0.001 0.000 A 0.00000 0.00000

For atomtype_FES.itp:
; atomtype_FES.atp
; Custom atom types for Fe and SF atoms
[ atomtypes ]
; name at.num mass charge ptype sigma (nm) epsilon (kJ/mol)
Fe 26 55.84500 0.6004 A 1.40000 0.01571
SF 16 32.06500 -0.5164 A 2.00000 0.25000

3

Topol of chain C and D are generated by pdb2gmx, but I replace the target Cys H bond with dummy atom.

The FES.itp:
[ moleculetype ]
; Name nrexcl
FES 3

[ atoms ]
; id type resnr residue atom cgnr charge mass
1 SF 1 FES S1 1 -0.51640 32.065
2 SF 1 FES S2 2 -0.51640 32.065
3 Fe 1 FES FE1 3 0.60040 55.845
4 Fe 1 FES FE2 4 0.60040 55.845

[ bonds ]
; i j func b0 (nm) kb (kJ/mol nm^2)
1 3 1 0.22 39080
2 4 1 0.22 39080
1 4 1 0.22 39080
2 3 1 0.22 39080

[ angles ]
; i j k func th0 (deg) cth (kJ/mol rad^2)
1 3 2 1 103 72.312
1 4 2 1 103 72.312
3 1 4 1 77 75.192
3 2 4 1 77 75.192

4

Finally is the CYS-Fes.itp that I want to define the interaction between my cluster to the Cys residue:
[ bonds ]
; i j func b0 (nm) kb (kJ/mol nm^2)
FE1 SG3 1 0.22 34559
FE1 SG4 1 0.22 34559
FE2 SG1 1 0.22 34559
FE2 SG2 1 0.22 34559

[ angles ]
; i j k func th0 (deg) cth (kJ/mol rad^2)
SG3 FE1 SG4 1 108 273.967
SG1 FE2 SG2 1 108 273.967
S1 FE1 SG3 1 103 72.312
S1 FE1 SG4 1 103 72.312
S2 FE2 SG1 1 103 72.312
S2 FE2 SG2 1 103 72.312

[ pairs ]
; i j func
FE1 SG3 1
FE1 SG4 1
FE2 SG1 1
FE2 SG2 1

But it usually return warning that: Too few parameter …

5

Now this system can be run by the software. But the result is not reasonable. The parameter are obtained frome papers.

6

Hey ,
If you have somehow figured it then would you mind sharing the parameters and itp file