Page 1 of 3

Print dipol vectors

Posted: Thu Oct 10, 2013 4:33 am
by michaelschauperl
Hi,

I just started using OpenMM and would need some help by setting up easy things.
Can you tell me how i can set up a reporter who is writing my the induced dipol moments(Amoeba) of an atom every x steps?

Thanks,

Michael

Re: Print dipol vectors

Posted: Thu Oct 10, 2013 2:09 pm
by peastman
Hi Michael,

Unfortunately, there's no API to query the induced dipole moments of individual atoms. It would be a reasonable thing to add, though. Enter it into the feature request tracker (under "Advanced" at the left side of this page) and I'll see if we can get it into the next version.

In the mean time, it's still possible if you don't mind hacking the source code. I can provide details if you're interested in trying that.

Peter

Re: Print dipol vectors

Posted: Fri Oct 11, 2013 1:32 am
by michaelschauperl
Hi Peter,

If you can provide some details on that for me, it would be great and i give it a go.

Thanks

Re: Print dipol vectors

Posted: Mon Oct 14, 2013 10:24 am
by peastman
Hi Michael,

To start with, if you haven't already done so, take a look at the developer guide. It gives information about the overall architecture and about how plugins (such as AMOEBA) are implemented.

As your starting point, take a look at AmoebaMultipoleForce::getSystemMultipoleMoments(). That computes the overall multipole moments of the whole system, which involves querying the individual induced dipoles, so you can see how that works. AmoebaMultipoleForce::getSystemMultipoleMoments() invokes getSystemMultipoleMoments() on its AmoebaMultipoleForceImpl, which in turn invokes getSystemMultipoleMoments() on its CalcAmoebaMultipoleForceKernel. So take a look in AmoebaCudaKernels.cpp and find CudaCalcAmoebaMultipoleForceKernel::getSystemMultipoleMoments().

The critical lines in that are line 1789, which ensures that the multipole moments have been calculated:

Code: Select all

ensureMultipolesValid(context);
and then line 1723, which downloads them from the GPU:

Code: Select all

inducedDipole->download(inducedDipoleVec);
inducedDipoleVec now contains all the atom dipole moments. The first three elements are the x, y, and z components for the first atom, the next three are the second atom, etc.

Does that make sense?

Peter

Re: Print dipol vectors

Posted: Thu Oct 17, 2013 3:29 pm
by peastman
Ok, I went ahead and implemented this myself. If you check out the latest code from github, you'll find that AmoebaMultipoleForce now has a method called getInducedDipoles(). Try it out and see if it works for you.

Peter

Re: Print dipol vectors

Posted: Fri Oct 18, 2013 12:51 am
by michaelschauperl
Hi Peter,

Thanks for including that feature. Unfortunately i was quite busy this week, so i wasn't able to try.
I will try it as soon as possible.
But a second question occured this week. Is it possible to add a charge plasma or a background charge grid to a OpenMM simulation?

Michael

Re: Print dipol vectors

Posted: Fri Oct 18, 2013 9:49 am
by peastman
Could you describe what you mean? If you're using PME, it already ignores the 0 frequency component, which is equivalent to filling all of space with a uniform charge density to neutralize the system. Is that what you're talking about, or something else?

Peter

Re: Print dipol vectors

Posted: Mon Oct 21, 2013 2:06 am
by michaelschauperl
That was exactly what i meant. Thanks for the quick response.

I just tried to do what you suggested to print out the dipol vectors, but i am not sure if i am doing everything in a sensible way.
Can you explain how i create a AmoebaMultipoleForce object in a sensible way. so what arguments i have to pass and so on. because if i create an AMF object and then call for example getNumMultipoles i get always 0)

Re: Print dipol vectors

Posted: Mon Oct 21, 2013 9:48 am
by peastman
I assume you're creating your system with ForceField.createSystem(), specifying the AMOEBA XML file? If so, then your system already contains an AmoebaMultipoleForce. You can get it like this:

Code: Select all

multipole = [f for f in system.getForces() if isinstance(f, AmoebaMultipoleForce)][0]
Peter

Re: Print dipol vectors

Posted: Thu Oct 24, 2013 4:03 am
by michaelschauperl
Thanks for your help so far. I am not sure what that third argument
std::vector<(Vec3,std::allocator<(Vec3)>)> dipoles
should be and how i define it.

Thanks for your help again.

Michael