Getting the charge of a single atom, per loop in MD Analysis

Question

I have been trying to use the partial charge of one particular ion to go through a calculation within mdanalysis. I have tried(This is just a snippet from the code that I know is throwing the error):

Cl = u.select_atoms('resname CLA and prop z <= 79.14') 
Lz = 79.14                          #Determined from system set-up
Q_sum = 0
COM = 38.42979431152344             #Determined from VMD
file_object1 = open(fors, 'a')
print(dcd, file = file_object1)
for ts in u.trajectory[200:]:
    frame = u.trajectory.frame
    time = u.trajectory.time
    for coord in Cl.positions:
        q= Cl.total_charge(Cl.position[coord][2])
        coords = coord - (Lz/COM)
        q_prof = q * (coords + (Lz / 2)) / Lz
        Q_sum = Q_sum + q_prof
        print(q)

But I keep getting an error associated with this.

How would I go about selecting this particular atom as it goes through the loop to get the charge of it in MD Analysis? Before I was setting q to equal a constant and the code ran fine so I know it is only this line that is throwing the error:

q = Cl.total_charge(Cl.position[coord][2])

Thanks for the help!

[total_charge()](https://docs.mdanalysis.org/stable/documentation_pages/core/groups.html?highlight=total_charge#MDAnalysis.core.groups.AtomGroup.total_charge) is a method of an AtomGroup and the only argument it takes (see docs) determines over which parts of the group the charge is summed. Can you descrive in words what you want to accomplish, then we can tell you _how_ to do it in MDAnalysis. — orbeckst, Aug 10 '21 at 15:54
Additionally, you can get the center of mass with the [center_of_mass()](https://docs.mdanalysis.org/stable/documentation_pages/core/groups.html?highlight=center_of_mass#MDAnalysis.core.groups.AtomGroup.center_of_mass) method of an AtomGroup. — orbeckst, Aug 10 '21 at 16:03

score 0 · Answer 1 · answered Dec 16 '21 at 03:23

I figured it out with:

def Q_code(dcd, topo):
Lz = u.dimensions[2]
Q_sum = 0
count = 0
CLAs = u.select_atoms('segid IONS or segid PROA or segid PROB or segid MEMB')
ini_frames = -200
n_frames = len(u.trajectory[ini_frames:])
for ts in u.trajectory[ini_frames:]:
    count += 1
    membrane = u.select_atoms('segid PROA or segid PROB or segid MEMB')
    COM = membrane.atoms.center_of_mass()[2]
    q_prof = CLAs.atoms.charges * (CLAs.positions[:,2] + (Lz/2 - COM))/Lz
    Q_instant = np.sum(q_prof)
    Q_sum += Q_instant
Q_av = Q_sum / n_frames
with open('Q_av.txt', 'a') as f:
    print('The Q_av for {}   is   {}'.format(s, Q_av), file = f)
return Q_av

Getting the charge of a single atom, per loop in MD Analysis

1 Answers1