Skip to content

Fix for systems with PDE and standard potential regions 馃挴

Peter Reinholdt requested to merge fix-stdembpot-and-pde into master

Description

Fixes missing polarizabilities in the dummy potential file. Polarizabilities may be present due to

  1. loprop regions
  2. standard_potential regions, which are polarizable (for example SEP).

The previous dummy potential writer:

def write_dummy_potential(system, filename):
    atom2site = {}
    site2atom = {}
    site_index = 1
    for region in system.regions.values():
        for fragment in region.fragments.values():
            for atom in fragment.atoms:
                site = Potential()
                system.potential[site_index] = site
                atom2site[atom.number] = site_index
                site2atom[site_index] = atom.number
                site.coordinate = atom.coordinate
                site.element = atom.element
                site_index += 1
    for region in system.regions.values():
        if region.use_polarizabilities:
            for fragment in region.fragments.values():
                for atom in fragment.atoms:
                    site = system.potential[atom2site[atom.number]]
                    # dummy polarizability
                    site.P11 = [1., 0., 0., 1., 0., 1.]
    system.write_potential(filename)

checked for polarizabilities by the region.use_polarizabilities property, which covers option 1). However, this property is not set for standard_potential regions, which led to dimension mismatches with PDE. To fix this, we can check if a region is either has region.use_polarizabilities or uses standard_potentials which are polarizable:

    ...
    for region in system.regions.values():
        if region.use_polarizabilities or\
        (region.use_standard_potentials and 'P11' in list(list(read_standard_potential(region.standard_potential_model).values())[0].values())[0].keys()):
            for fragment in region.fragments.values():
                for atom in fragment.atoms:
                    site = system.potential[atom2site[atom.number]]
                    # dummy polarizability
                    site.P11 = [1., 0., 0., 1., 0., 1.]
    system.write_potential(filename)

Types of changes

  • Bug fix (non-breaking change which fixes an issue)
Edited by Peter Reinholdt

Merge request reports

Loading