There was a really interesting question on SO (http://stackoverflow.com/questions/15709496) that grabbed my attention. Unfortunately, it has been closed, but I thought it was a great question, and had started a mash-up answer. Seeing as it is now closed, I will work through my model here.
We have an input file that contains the measured (x, y, z) coordinates of all three atoms in a water molecule. Our job is to extend the (x, y, z) data to represent any number of water molecules. Coordinates should be generated psuedo-randomly and no molecule should be close that 0.2 nm from any other.
This project is open to anyone who wants to help out! Feel free to post comments, suggestions or different implementations.
Below is my initial work on this (it is a work in progress):
These functions do the following:
- import_data(infile): reads the measured data from the file and formats it
- data_fields(data): if necessary, this function will group all ‘like’ fields (i.e. create a list of all x coords, or all atom_IDs)(and below…)
- distance(a, b): calculates the distance between two 3-dimensional points
- gen_coords(): this one is interesting, it takes as input the coordinates of a neighbouring atom, and then adds a psuedo-random shift to create the new coords. The size of the shift is determined by (min_range, max_range) which is approx. the distance of a covalent or polar bond as necessary
- test_coords() test whether the psuedo-random shift above is within an appropriate range (i.e. covalent bonds should be approx. 1 Armstrong in length, and polar bonds are approximately 2 Armstrongs long). Will recall gen_coords() until condition is satisfied
Still to do:
- generate coords for HW2. The above works well for O and HW1, but is not accurate for HW2. To be more accurate, HW2 should be a 120 degree rotation of HW1. Ideally, this would be accomplished by a rotation matrix.
- Add molecule coords and meta-data back into the datafile.
- Lastly, plot it, see what the model looks like.
Writing out to file and plotting are straightforward. The rotation matrix is a little trickier, and will require me to re-learn some maths I haven’t used for years.
This project is open to anyone who would like to help out – feel free to post comments, suggestions or implementations on this 🙂