Periodic bond rendering

Development
1

On Sun, Apr 15, 2012 at 11:22 AM, Tuukka Verho tuukka.verho@aalto.fi wrote:

Your screenshot looks quite nice, indeed multiple intersections are correctly
handled. Though I’m getting the feeling that the code is quite large, and the
support for coloring and higher order bonds still need to be added as you say.
While thinking of a more straightforward way to achieve the same result I
stumbled upon some example code to do boolean operations with opengl

http://www.hackchina.com/en/r/78488/csg.c__html

If the clipping could be done with a logical AND operation, multiple bonds
etc. would be handled without any extra effort. Also cases when the cylinder
intersects with not just one but two or three walls would be nicely handled.

I will experiment a bit and let you know if anything comes out of it.

I would advise against investing too much effort in this solution.
While it certainly should work for OpenGL, all rendering should be
done through painters so that alternative rendering backends (POV-ray,
etc) can have the same features. Clipping the primitives through GL
won’t have an effect on the other backends, which is why I’ve added
the extra code to actually draw the clipped primitives in the painter.

Dave

2

On Sunday, April 15, 2012 11:47:25 AM David Lonie wrote:

On Sun, Apr 15, 2012 at 11:22 AM, Tuukka Verho tuukka.verho@aalto.fi wrote:

Your screenshot looks quite nice, indeed multiple intersections are
correctly handled. Though I’m getting the feeling that the code is quite
large, and the support for coloring and higher order bonds still need to
be added as you say. While thinking of a more straightforward way to
achieve the same result I stumbled upon some example code to do boolean
operations with opengl

http://www.hackchina.com/en/r/78488/csg.c__html

If the clipping could be done with a logical AND operation, multiple bonds
etc. would be handled without any extra effort. Also cases when the
cylinder intersects with not just one but two or three walls would be
nicely handled.

I will experiment a bit and let you know if anything comes out of it.

I would advise against investing too much effort in this solution.
While it certainly should work for OpenGL, all rendering should be
done through painters so that alternative rendering backends (POV-ray,
etc) can have the same features. Clipping the primitives through GL
won’t have an effect on the other backends, which is why I’ve added
the extra code to actually draw the clipped primitives in the painter.

Dave

Ok let’s stay with a general solution.

  • Tuukka
3

Cool – feel free to take a look at what I’ve got and let me know if
you can spot any possible improvements:

https://github.com/dlonie/avogadro/tree/ENH_intercell_bonds

Hi Dave,

Your algorithm for finding the shortest bond vector doesn’t work right for unit
cells with angles other than 90deg. I changed it to something that should
yield correct results but I can’t compile it because I get:

/home/tuukka/code/dlonie/libavogadro/src/engines/bsdyengine.cpp:35:34: fatal
error: avogadro/obeigenconv.h: No such file or directory

However I made a patch of my changes, you can see if it compiles.

  • Tuukka
4

Hi Tuukka,

On Sun, Apr 15, 2012 at 1:58 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

Hi Dave,

Your algorithm for finding the shortest bond vector doesn’t work right for unit
cells with angles other than 90deg.

You’re right! I saw that method in a patch that went into openbabel
and thought it was clever – mapping the end position into a
Wigner-Seitz cell around the start position is far less expensive than
checking 27 images explicitly. But yes, this will only work for cubic
lattices. Thanks for catching that – I’ll patch up OpenBabel with a
more reliable method, too.

I changed it to something that should
yield correct results but I can’t compile it because I get:

/home/tuukka/code/dlonie/libavogadro/src/engines/bsdyengine.cpp:35:34: fatal
error: avogadro/obeigenconv.h: No such file or directory

Oops – This branch does depends on this old patch on gerrit:

http://review.source.kitware.com/#/c/4238/

Maybe one of the other core devs will find some time to review the
backlog of patches soon? For now, just copying the obeigenconv.h
header from libavogadro/src/extensions/crystallography to some place
that bsdyengine.cpp can find it will let the code compile. I had a
stale copy in my working tree, so I hadn’t noticed it missing from the
patch :wink:

However I made a patch of my changes, you can see if it compiles.

Thanks – I tweaked it a little to make it a bit more efficient and
use Eigen vectors/matrices instead of OB’s linear algebra
implementation. I’ve put the new version up on github.

Dave

5

Your algorithm for finding the shortest bond vector doesn’t work right for
unit cells with angles other than 90deg.

You’re right! I saw that method in a patch that went into openbabel
and thought it was clever – mapping the end position into a
Wigner-Seitz cell around the start position is far less expensive than
checking 27 images explicitly. But yes, this will only work for cubic
lattices. Thanks for catching that – I’ll patch up OpenBabel with a
more reliable method, too.

At some point I also felt tempted to just do fractionalVector[i] %= .5 but I
got the feeling that it’s not very rigorous :)…

However I made a patch of my changes, you can see if it compiles.

Thanks – I tweaked it a little to make it a bit more efficient and
use Eigen vectors/matrices instead of OB’s linear algebra
implementation. I’ve put the new version up on github.

I tidied it up a bit more and changed some of the variable names to better
ones. I think this version is also more efficient if you count the function
calls. I managed to compile this time, thanks…

Drawing the bonds in several pieces may lead to funny effects in some cases,
there’s an example in the attachments. I don’t know how big a problem that is
but some such cases might end up being reported as bugs…

  • Tuukka
6

Hi Tuukka,

On Mon, Apr 16, 2012 at 9:59 AM, Tuukka Verho tuukka.verho@aalto.fi wrote:

I tidied it up a bit more and changed some of the variable names to better
ones. I think this version is also more efficient if you count the function
calls. I managed to compile this time, thanks…

The original code was more efficient. Pretty much all of the Eigen
code (MatrixBase::col(), operator*/+, etc) end up getting inlined by
the compiler, so there is no extra overhead from traversing the stack.
I haven’t checked the assembly, but I believe it ends up being
equivalent to just working with straight arrays. By default, Eigen
also stores matrices as column-major, so operations on
MatrixBase.col() should be able to be vectorized, as well.

The cellMatrixT * q takes 9 multiplications and is performed for each
of the 27 image vectors, so this costs 243 FLOPs per bond. Contrasted
with caching the vectors for each loop, the cost is reduced to 3 * 3
(a) + 9 * 3 (b) + 27 * 3 (c) = 117 FLOPs per bond, less than half of
the required operations for the matrix approach. The number of cheaper
additions should be roughly the same. I’m ignoring the squared norm
calculation (3 FLOPs), but these are unavoidable.

That said, my experience in performance analysis is limited, and I’m
not sure how the caching overhead may affect things. I assume the
cached vectors can be stored in registers, or L1 cache at worst, so
it’s probably rather trivial.

Drawing the bonds in several pieces may lead to funny effects in some cases,
there’s an example in the attachments. I don’t know how big a problem that is
but some such cases might end up being reported as bugs…

Those are nasty looking indeed! We should try to think of a way to
clean those up before anything gets submitted. The main problem here
is that the clippedCylinder only clips along one plane, so if a bond
segment is near two planes part of the cylinder will jut out. A
potential solution would be to have the clippedCylinder take a set of
planes as arguments, but this would be a pain to code.

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

Dave

7

David Lonie kirjoitti 16.4.2012 kello 18.39:

The original code was more efficient. Pretty much all of the Eigen
code (MatrixBase::col(), operator*/+, etc) end up getting inlined by
the compiler, so there is no extra overhead from traversing the stack.
I haven’t checked the assembly, but I believe it ends up being
equivalent to just working with straight arrays. By default, Eigen
also stores matrices as column-major, so operations on
MatrixBase.col() should be able to be vectorized, as well.

Ok. But please change the variable names, the original ones were
confusing. And there’s no need to calculate endPos-startPos again each
time.

Drawing the bonds in several pieces may lead to funny effects in
some cases,
there’s an example in the attachments. I don’t know how big a
problem that is
but some such cases might end up being reported as bugs…

Those are nasty looking indeed! We should try to think of a way to
clean those up before anything gets submitted. The main problem here
is that the clippedCylinder only clips along one plane, so if a bond
segment is near two planes part of the cylinder will jut out. A
potential solution would be to have the clippedCylinder take a set of
planes as arguments, but this would be a pain to code.

With up to 3 clipping planes there are so many possibilities and
corner cases that it is virtually impossible to avoid all artifacts.
We just have to settle how many cases we want to cover (we don’t want
to end up in a situation where 90% of drawing code in Avogadro is
related to this). I was actually a bit unsure whether drawing those
“isolated” bond segments is a good idea because they often end up
being clipped in complicated ways.

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

I would vote for this option…

  • Tuukka
8

On Mon, Apr 16, 2012 at 12:50 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

David Lonie kirjoitti 16.4.2012 kello 18.39:

With up to 3 clipping planes there are so many possibilities and corner
cases that it is virtually impossible to avoid all artifacts. We just have
to settle how many cases we want to cover (we don’t want to end up in a
situation where 90% of drawing code in Avogadro is related to this). I was
actually a bit unsure whether drawing those “isolated” bond segments is a
good idea because they often end up being clipped in complicated ways.

Is this going into a specialized plugin, or straight into bsdyengine?
I didn’t look at the code, but it seems pretty specialized and I
wondered what (if any) impact this would have where there is no
periodic boundary.

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

I would vote for this option…

I think this would be reasonable, and I think virtually all rendering
engines have some form of clipping support that we could map to.

Marcus

9

On Mon, Apr 16, 2012 at 12:50 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

Ok. But please change the variable names, the original ones were confusing.

Sure – which ones are confusing? They seem straightforward to me:

[abc]Trans: Translation vector along the a, b, or c axis
curImage: Vector pointing to the current translational image of the atoms/bonds

What would you suggest instead?

And there’s no need to calculate endPos-startPos again each time.

Good point, I’ll make that change.

With up to 3 clipping planes there are so many possibilities and corner
cases that it is virtually impossible to avoid all artifacts. We just have
to settle how many cases we want to cover (we don’t want to end up in a
situation where 90% of drawing code in Avogadro is related to this). I was
actually a bit unsure whether drawing those “isolated” bond segments is a
good idea because they often end up being clipped in complicated ways.

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

I would vote for this option…

I think for now we should add an option to the bsdyengine settings
widget that will allow users to switch off isolated segment rendering,
or periodic bonding entirely. I can work on this later if no one gets
to it before I have time.

Dave

10

On Mon, Apr 16, 2012 at 1:13 PM, David Lonie loniedavid@gmail.com wrote:

On Mon, Apr 16, 2012 at 12:50 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

Ok. But please change the variable names, the original ones were confusing.

Sure – which ones are confusing? They seem straightforward to me:

[abc]Trans: Translation vector along the a, b, or c axis
curImage: Vector pointing to the current translational image of the atoms/bonds

What would you suggest instead?

And there’s no need to calculate endPos-startPos again each time.

Good point, I’ll make that change.

With up to 3 clipping planes there are so many possibilities and corner
cases that it is virtually impossible to avoid all artifacts. We just have
to settle how many cases we want to cover (we don’t want to end up in a
situation where 90% of drawing code in Avogadro is related to this). I was
actually a bit unsure whether drawing those “isolated” bond segments is a
good idea because they often end up being clipped in complicated ways.

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

I would vote for this option…

I think for now we should add an option to the bsdyengine settings
widget that will allow users to switch off isolated segment rendering,
or periodic bonding entirely. I can work on this later if no one gets
to it before I have time.

It would be nice to contain a lot of these changes in a more
specialized plugin, but if they are at least well isolated that should
be OK.

Marcus

11

On Mon, Apr 16, 2012 at 12:56 PM, Marcus D. Hanwell
mhanwell@gmail.com wrote:

Is this going into a specialized plugin, or straight into bsdyengine?
I didn’t look at the code, but it seems pretty specialized and I
wondered what (if any) impact this would have where there is no
periodic boundary.

I’ll save you the trouble of digging through the patch :wink:

https://github.com/dlonie/avogadro/blob/ENH_intercell_bonds/libavogadro/src/engines/bsdyengine.cpp#L150

If there’s no unit cell defined, the bonds are rendered as usual.
Otherwise each bond is shortened if needed, and tested to see if it
intersects any face of the plane.

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

I would vote for this option…

I think this would be reasonable, and I think virtually all rendering
engines have some form of clipping support that we could map to.

Great, should we slot this in for 2.0, or should we go ahead and add
it now? I think we could do it without any API breaks.

My only concern with this approach is that it won’t close off the
cylinders, leading to odd backface culling artifacts.

Dave

12

On Mon, Apr 16, 2012 at 1:18 PM, Marcus D. Hanwell
mhanwell@gmail.com wrote:

It would be nice to contain a lot of these changes in a more
specialized plugin, but if they are at least well isolated that should
be OK.

I was trying to think of a way to implement this in a maintainable way
across all of the engines, maybe with a MoleculeRepresentation (or
similar) layer between Molecule and Engine, but a new engine may be
the best way to go for now. Or add options to the engine settings so
that users could select how to handle periodic bonds:

  • Expanded: Complete all molecules that exceed the unit cell:
    (see the graphical abstract on http://pubs.acs.org/doi/abs/10.1021/ja201786y )
  • Clipped: Think of a 3x3x3 supercell clipped so that only fragments
    inside the central cell are shown:
    Zink Crystal Video It’s a video, but the preview frame
    shows what I mean. This is what I’m trying to get implemented at the
    moment.
  • Simple periodic: Like Tuukka’s original patch, which just renders
    the bond halves, not clipping them at the boundaries.

I prefer adding options to the existing engines. It’s just a shame to
have to recalculate all of the translated geometry on every render for
each engine, and it will suck to implement and maintain this across
all engines :-/

Dave

13

David Lonie kirjoitti 16.4.2012 kello 20.13:

On Mon, Apr 16, 2012 at 12:50 PM, Tuukka Verho
tuukka.verho@aalto.fi wrote:

Ok. But please change the variable names, the original ones were
confusing.

Sure – which ones are confusing? They seem straightforward to me:

[abc]Trans: Translation vector along the a, b, or c axis
curImage: Vector pointing to the current translational image of the
atoms/bonds

What would you suggest instead?

No I meant the ones I changed in my patch:

sqNormDifference → minLengthSq
curSqDist → lengthSq

14

On Mon, Apr 16, 2012 at 2:11 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

No I meant the ones I changed in my patch:

sqNormDifference → minLengthSq
curSqDist → lengthSq

Ok, done.

Dave

15

David Lonie kirjoitti 16.4.2012 kello 20.19:

Maybe Avogadro 2.0 painters could have a set of clipping planes so
that we could more easily draw these structures? Then we could just
use native GL clipping and emulate it somehow for the other
paintdevices…

I would vote for this option…

I think this would be reasonable, and I think virtually all rendering
engines have some form of clipping support that we could map to.

Great, should we slot this in for 2.0, or should we go ahead and add
it now? I think we could do it without any API breaks.

My only concern with this approach is that it won’t close off the
cylinders, leading to odd backface culling artifacts.

Dave

No, you can also close the cylinders. First you draw the pixels of the
cylinder that are inside the clipping box, then you draw the parts of
the clipping box that are inside the cylinder. Voila!

  • Tuukka
16

On Mon, Apr 16, 2012 at 2:42 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

David Lonie kirjoitti 16.4.2012 kello 20.19:

My only concern with this approach is that it won’t close off the
cylinders, leading to odd backface culling artifacts.

No, you can also close the cylinders. First you draw the pixels of the
cylinder that are inside the clipping box, then you draw the parts of the
clipping box that are inside the cylinder. Voila!

Ah, good point :slight_smile: That will work for the GL backend. Not sure if the
other renderers will support non-planar clip operations, but we can
always emulate it if needed.

Dave

17

David Lonie kirjoitti 16.4.2012 kello 20.37:

On Mon, Apr 16, 2012 at 1:18 PM, Marcus D. Hanwell
mhanwell@gmail.com wrote:

It would be nice to contain a lot of these changes in a more
specialized plugin, but if they are at least well isolated that
should
be OK.

I was trying to think of a way to implement this in a maintainable way
across all of the engines, maybe with a MoleculeRepresentation (or
similar) layer between Molecule and Engine, but a new engine may be
the best way to go for now.

Maybe the clipping code could be simply taken out from the engines? We
could have a more sophisticated drawCylinder() method that would do
the clipping etc. and also draw two-colored cylinders? On the other
hand then it would need to be replicated in each painter backend
separately – unless the backend independent stuff would be in the
Painter superclass for example.

Or add options to the engine settings so
that users could select how to handle periodic bonds:

This is certainly an interesting option (though you need to be careful
with infinite molecules such as covalent crystals). This requires a
completely different approach – maybe Molecule could have a method to
generate all the atoms and bonds needed to create one complete
instance of itself (in the sense that all parts in the primary cell
are expanded).

  • Clipped: Think of a 3x3x3 supercell clipped so that only fragments
    inside the central cell are shown:
    Zink Crystal Video It’s a video, but the preview frame
    shows what I mean. This is what I’m trying to get implemented at the
    moment.
  • Simple periodic: Like Tuukka’s original patch, which just renders
    the bond halves, not clipping them at the boundaries.

I prefer adding options to the existing engines. It’s just a shame to
have to recalculate all of the translated geometry on every render for
each engine, and it will suck to implement and maintain this across
all engines :-/

Wait a minute… I’m actually new to opengl, but do you really repaint
the scene on each render? I mean can’t you store the opengl commands
and execute them just with a different transformation matrix when the
view is e.g. rotated?

  • Tuukka
18

Hi Tuukka,

On Mon, Apr 16, 2012 at 2:56 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

David Lonie kirjoitti 16.4.2012 kello 20.37:

I was trying to think of a way to implement this in a maintainable way
across all of the engines, maybe with a MoleculeRepresentation (or
similar) layer between Molecule and Engine, but a new engine may be
the best way to go for now.

Maybe the clipping code could be simply taken out from the engines? We could
have a more sophisticated drawCylinder() method that would do the clipping
etc. and also draw two-colored cylinders? On the other hand then it would
need to be replicated in each painter backend separately – unless the
backend independent stuff would be in the Painter superclass for example.

This is similar to the idea of putting clipping plane support directly
into the painters, which I think we’ve agreed would be a good move.
This would allow us to clip spheres, lines, etc, too for atoms and
wireframe in addition to just the cylinders.

Or add options to the engine settings so
that users could select how to handle periodic bonds:

This is certainly an interesting option (though you need to be careful with
infinite molecules such as covalent crystals). This requires a completely
different approach – maybe Molecule could have a method to generate all the
atoms and bonds needed to create one complete instance of itself (in the
sense that all parts in the primary cell are expanded).

I’ve got some code in XtalOpt that does something similar when the
time comes. Nothing fancy, it just walks through the bonds, shortening
them as it goes. I’ve haven’t tried it on infinite systems, but it
shouldn’t be too hard to adapt.

  • Clipped: Think of a 3x3x3 supercell clipped so that only fragments
    inside the central cell are shown:
    Zink Crystal Video It’s a video, but the preview frame
    shows what I mean. This is what I’m trying to get implemented at the
    moment.
  • Simple periodic: Like Tuukka’s original patch, which just renders
    the bond halves, not clipping them at the boundaries.

I prefer adding options to the existing engines. It’s just a shame to
have to recalculate all of the translated geometry on every render for
each engine, and it will suck to implement and maintain this across
all engines :-/

Wait a minute… I’m actually new to opengl, but do you really repaint the
scene on each render? I mean can’t you store the opengl commands and execute
them just with a different transformation matrix when the view is e.g.
rotated?

Yes, the renderOpaque and renderTransparent methods are run each
render (which is why I was so concerned about efficiency earlier :wink:
). AIUI, this is something that will be improved for Avogadro 2.0.

Dave

19

This is similar to the idea of putting clipping plane support directly
into the painters, which I think we’ve agreed would be a good move.
This would allow us to clip spheres, lines, etc, too for atoms and
wireframe in addition to just the cylinders.

Indeed, many people like to clip an orbital or other surface.

-Geoff

20

Hi Dave,

On Tuesday, April 17, 2012 01:17:45 PM David Lonie wrote:

Hi Tuukka,

On Mon, Apr 16, 2012 at 2:56 PM, Tuukka Verho tuukka.verho@aalto.fi wrote:

David Lonie kirjoitti 16.4.2012 kello 20.37:

I was trying to think of a way to implement this in a maintainable way
across all of the engines, maybe with a MoleculeRepresentation (or
similar) layer between Molecule and Engine, but a new engine may be
the best way to go for now.

Maybe the clipping code could be simply taken out from the engines? We
could have a more sophisticated drawCylinder() method that would do the
clipping etc. and also draw two-colored cylinders? On the other hand then
it would need to be replicated in each painter backend separately –
unless the backend independent stuff would be in the Painter superclass
for example.
This is similar to the idea of putting clipping plane support directly
into the painters, which I think we’ve agreed would be a good move.
This would allow us to clip spheres, lines, etc, too for atoms and
wireframe in addition to just the cylinders.

That’s right… But I was thinking a bit further: the painter could also take
care of calculating the shortest distance, multiple intersections and so on.
So if the engine wants to wrap the cylinders it just passes the unit cell to
the painter and everything just works. In that scheme the current way of
drawing a two-colored bond with two separate drawCylinder calls wouldn’t work,
so there should be a drawTwoColoredCylinder method in the painter.

  • Tuukka