2
$\begingroup$

I am trying to build a supercell for molecular crystals from a .cif file. With the ASE package, the supercell has the perfect PBC, but the atom sequence is very random and might be difficult to match with the force field topology to run the simulation in GROMACS. With Mercury, it extends the cell but doesn't provide the cell size for PBC. If I multiply the cell size just by the original cell size, it does not fulfill the PBC condition. Do you have any suggestions or other tools to fix this? I appreciate your help.

####################################################################### 
# 
# This file contains crystal structure data downloaded from the 
# Cambridge Structural Database (CSD) hosted by the Cambridge 
# Crystallographic Data Centre (CCDC).
# 
# Full information about CCDC data access policies and citation 
# guidelines are available at http://www.ccdc.cam.ac.uk/access/V1 
# 
# Audit and citation data items may have been added by the CCDC. 
# Please retain this information to preserve the provenance of 
# this file and to allow appropriate attribution of the data. 
# 
#######################################################################

data_urea_310
_audit_block_doi                 10.5517/ccdc.csd.ccsknnm
_database_code_depnum_ccdc_archive 'CCDC 731961'
loop_
_citation_id
_citation_doi
_citation_year
1 10.1021/jp904942c 2009
_audit_update_record             
;
2009-05-13 deposited with the CCDC. 2023-11-08 downloaded from the CCDC.
;

_audit_creation_method           SHELXL-97
_chemical_name_systematic        
;
urea
;
_chemical_name_common            urea
_chemical_melting_point          405
_chemical_formula_moiety         'C H4 N2 O'
_chemical_formula_sum            'C H4 N2 O'
_chemical_formula_weight         60.06

_symmetry_cell_setting           orthorhombic
_symmetry_space_group_name_H-M   'P 21 21 2'
_symmetry_space_group_name_Hall  'P 2 2ab'

loop_
_symmetry_equiv_pos_as_xyz
'x, y, z'
'-x, -y, z'
'-x+1/2, y+1/2, -z'
'x+1/2, -y+1/2, -z'

_cell_length_a                   3.414(3)
_cell_length_b                   7.360(8)
_cell_length_c                   4.606(10)
_cell_angle_alpha                90.00
_cell_angle_beta                 90.00
_cell_angle_gamma                90.00
_cell_volume                     115.7(3)
_cell_formula_units_Z            2
_cell_measurement_temperature    296(2)
_cell_measurement_pressure       3100000
_cell_measurement_reflns_used    95
_cell_measurement_theta_min      5.54
_cell_measurement_theta_max      26.59

_exptl_crystal_description       plate
_exptl_crystal_colour            colourless
_exptl_crystal_size_max          0.37
_exptl_crystal_size_mid          0.37
_exptl_crystal_size_min          0.10
_exptl_crystal_density_meas      ?
_exptl_crystal_density_diffrn    1.723
_exptl_crystal_density_method    'not measured'
_exptl_crystal_F_000             64
_exptl_absorpt_coefficient_mu    0.148
_exptl_absorpt_correction_type   integration
_exptl_absorpt_correction_T_min  0.69
_exptl_absorpt_correction_T_max  0.95
_exptl_absorpt_process_details   
;
Katrusiak, A. (2003). REDSHABS - Program for correcting
reflections intensities for DAC absorption, gasket shadowing
and sample crystal absorption. Adam Mickiewicz University, Pozna\'n.
Katrusiak, A. (2004). Z. Kristallogr. 219, 461-467
;

_exptl_special_details           
;
Data were collected at room temperature and pressure of 3.10(5) GPa
(3100000 kPa) with the crystal obtained by the in-situ high-pressure
crystallization technique. Pressure was determined by monitoring
the shift of the ruby R1-fluorescence line.
;

_diffrn_ambient_temperature      296(2)
_diffrn_radiation_wavelength     0.71073
_diffrn_radiation_type           MoK\a
_diffrn_radiation_source         'fine-focus sealed tube'
_diffrn_radiation_monochromator  graphite
_diffrn_measurement_device_type  'Kuma KM4CCD \k geometry'
_diffrn_measurement_method       
;HP omega scans - for more details see:
A. Budzianowski, A. Katrusiak in High-Pressure Crystallography
(Eds.: A. Katrusiak, P. F. McMillan),
Dordrecht: Kluwer Acad. Publ., 2004 pp.157-168
;
_diffrn_detector_area_resol_mean 16.4
_diffrn_reflns_number            327
_diffrn_reflns_av_R_equivalents  0.1123
_diffrn_reflns_av_sigmaI/netI    0.0705
_diffrn_reflns_limit_h_min       -4
_diffrn_reflns_limit_h_max       4
_diffrn_reflns_limit_k_min       -7
_diffrn_reflns_limit_k_max       8
_diffrn_reflns_limit_l_min       -2
_diffrn_reflns_limit_l_max       2
_diffrn_reflns_theta_min         5.54
_diffrn_reflns_theta_max         26.59
_diffrn_measured_fraction_theta_max 0.419
_diffrn_reflns_theta_full        26.59
_diffrn_measured_fraction_theta_full 0.419
_diffrn_standards_number         ?
_diffrn_standards_interval_count ?
_diffrn_standards_interval_time  ?
_diffrn_standards_decay_%        ?
_refine_special_details          
;
Refinement of F^2^ against ALL reflections. The weighted R-factor wR and
goodness of fit S are based on F^2^, conventional R-factors R are based
on F, with F set to zero for negative F^2^. The threshold expression of
F^2^ > 2sigma(F^2^) is used only for calculating R-factors(gt) etc. and is
not relevant to the choice of reflections for refinement. R-factors based
on F^2^ are statistically about twice as large as those based on F, and R-
factors based on ALL data will be even larger.
;
_reflns_number_total             95
_reflns_number_gt                69
_reflns_threshold_expression     >2sigma(I)
_refine_ls_structure_factor_coef Fsqd
_refine_ls_matrix_type           full
_refine_ls_number_reflns         95
_refine_ls_number_parameters     20
_refine_ls_number_restraints     0
_refine_ls_R_factor_all          0.0573
_refine_ls_R_factor_gt           0.0332
_refine_ls_wR_factor_ref         0.0839
_refine_ls_wR_factor_gt          0.0661
_refine_ls_goodness_of_fit_ref   1.217
_refine_ls_restrained_S_all      1.217
_refine_ls_shift/su_max          0.000
_refine_ls_shift/su_mean         0.000
_refine_ls_weighting_scheme      calc
_refine_ls_weighting_details     
'calc w=1/[\s^2^(Fo^2^)+(0.0244P)^2^+0.0087P] where P=(Fo^2^+2Fc^2^)/3'
_atom_sites_solution_primary     direct
_atom_sites_solution_secondary   difmap
_atom_sites_solution_hydrogens   geom
_refine_ls_hydrogen_treatment    constr
_refine_ls_extinction_method     none
_refine_ls_extinction_coef       ?
_refine_diff_density_max         0.206
_refine_diff_density_min         -0.196
_refine_diff_density_rms         0.059
_refine_ls_abs_structure_details 'Flack H D (1983), Acta Cryst. A39, 876-881'
_refine_ls_abs_structure_Flack   7(8)
loop_
_atom_type_symbol
_atom_type_description
_atom_type_scat_dispersion_real
_atom_type_scat_dispersion_imag
_atom_type_scat_source
C C 0.0033 0.0016 'International Tables Vol C Tables 4.2.6.8 and 6.1.1.4'
H H 0.0000 0.0000 'International Tables Vol C Tables 4.2.6.8 and 6.1.1.4'
N N 0.0061 0.0033 'International Tables Vol C Tables 4.2.6.8 and 6.1.1.4'
O O 0.0106 0.0060 'International Tables Vol C Tables 4.2.6.8 and 6.1.1.4'

_computing_data_collection       'CrysAlis (Oxford Diffraction, 2002)'
_computing_cell_refinement       'CrysAlis (Oxford Diffraction, 2002)'
_computing_data_reduction        'CrysAlis (Oxford Diffraction, 2002)'
_computing_structure_solution    'SHELXS-97 (Sheldrick, 1990)'
_computing_structure_refinement  'SHELXL-97 (Sheldrick, 1997)'
_computing_molecular_graphics    
;X-Seed (Barbour, 2001)
and POV-Ray (Persistence of Vision, 2004)
;
_computing_publication_material  'SHELXL-97 (Sheldrick, 1997)'

loop_
_atom_site_label
_atom_site_type_symbol
_atom_site_fract_x
_atom_site_fract_y
_atom_site_fract_z
_atom_site_U_iso_or_equiv
_atom_site_adp_type
_atom_site_occupancy
_atom_site_symmetry_multiplicity
_atom_site_calc_flag
_atom_site_refinement_flags
_atom_site_disorder_assembly
_atom_site_disorder_group
C1 C 0.0000 0.5000 0.870(4) 0.023(8) Uani 1 2 d S . .
O1 O 0.0000 0.5000 1.136(3) 0.046(6) Uani 1 2 d S . .
N1 N 0.1270(9) 0.6429(4) 0.7181(19) 0.044(2) Uani 1 1 d . . .
H1A H 0.2115 0.7379 0.8070 0.053 Uiso 1 1 calc R . .
H1B H 0.1244 0.6399 0.5315 0.053 Uiso 1 1 calc R . .

loop_
_atom_site_aniso_label
_atom_site_aniso_U_11
_atom_site_aniso_U_22
_atom_site_aniso_U_33
_atom_site_aniso_U_23
_atom_site_aniso_U_13
_atom_site_aniso_U_12
C1 0.023(3) 0.032(5) 0.01(2) 0.000 0.000 0.007(3)
O1 0.053(3) 0.036(3) 0.05(2) 0.000 0.000 0.001(2)
N1 0.0442(19) 0.0171(18) 0.071(8) 0.007(3) -0.003(3) -0.0087(18)

_geom_special_details            
;
All esds (except the esd in the dihedral angle between two l.s. planes)
are estimated using the full covariance matrix. The cell esds are taken
into account individually in the estimation of esds in distances, angles
and torsion angles; correlations between esds in cell parameters are only
used when they are defined by crystal symmetry. An approximate (isotropic)
treatment of cell esds is used for estimating esds involving l.s. planes.
;

loop_
_geom_bond_atom_site_label_1
_geom_bond_atom_site_label_2
_geom_bond_distance
_geom_bond_site_symmetry_2
_geom_bond_publ_flag
C1 O1 1.225(13) . ?
C1 N1 1.336(9) 2_565 ?
C1 N1 1.336(9) . ?
N1 H1A 0.8600 . ?
N1 H1B 0.8600 . ?

loop_
_geom_angle_atom_site_label_1
_geom_angle_atom_site_label_2
_geom_angle_atom_site_label_3
_geom_angle
_geom_angle_site_symmetry_1
_geom_angle_site_symmetry_3
_geom_angle_publ_flag
O1 C1 N1 121.6(6) . 2_565 ?
O1 C1 N1 121.6(6) . . ?
N1 C1 N1 116.8(12) 2_565 . ?
C1 N1 H1A 120.0 . . ?
C1 N1 H1B 120.0 . . ?
H1A N1 H1B 120.0 . . ?


#END
$\endgroup$
4
  • 1
    $\begingroup$ Could you share your cif file, I might try to help you out. $\endgroup$ Jan 9 at 2:52
  • 2
    $\begingroup$ I have added the cif file. if you want to have a try. Appreciate your help. $\endgroup$
    – Pradip Si
    Jan 10 at 16:35
  • $\begingroup$ Save coordinates of the unit cell in pdb, delete atoms that will overlap with others when replicated, use gromacs' genconf to replicate the unit cell. $\endgroup$
    – Anon
    Jan 10 at 20:15
  • 1
    $\begingroup$ Gromacs "genconf' does not replicate the cell properly. You can also find here why it can be problematic gcm.upc.edu/en/members/luis-carlos/molecular-dynamics/… $\endgroup$
    – Pradip Si
    Jan 10 at 21:04

1 Answer 1

2
$\begingroup$

I tried to produce a supercell using ASE while ensuring PBC to be True, and I fixed the atoms in the original unit cell (maybe this would fix the problem with having random atom sequence - I am not very sure about this part)

from ase import io, Atoms
from ase.constraints import FixAtoms
import numpy as np
structure = io.read('YourStructure.cif')
supercell_size = [2, 2, 2]
lattice_vectors = structure.cell
atomic_positions = structure.positions
supercell_atoms = []
for i in range(supercell_size[0]):
    for j in range(supercell_size[1]):
        for k in range(supercell_size[2]):
            supercell_atoms.extend(
                [atom for atom in atomic_positions + np.dot([i, j, k], lattice_vectors)]
            )
supercell = Atoms(symbols=structure.get_chemical_symbols() * np.prod(supercell_size),
                  positions=supercell_atoms,
                  cell=np.dot(np.diag(supercell_size), lattice_vectors),
                  pbc=True)
constraint = FixAtoms(indices=range(len(atomic_positions)))
supercell.set_constraint(constraint)
io.write('supercell.cif', supercell, format='cif')

But I have to mention that ASE (cif.py) was throwing the following warning:

crystal system 'orthorhombic' is not interpreted for space group Spacegroup(18, setting=1). This may result in wrong setting! 

So I am not so sure about the validation of the cif file you have.

$\endgroup$

You must log in to answer this question.