Re: [boost] interest on a geometry library

From: Janek Kozicki (janek_listy_at_[hidden])
Date: 2006-03-20 15:35:53

Andy Little said:     (by the date of Mon, 20 Mar 2006 14:14:28 -0000)

> > I will try to help you with that geometry library, because I want to
> > use boost too, instead of dragging everywhere similar library
> > (with quaternions, and all that 3D stuff). I'm mainly using that library
> > for discrete simulation modelling (shameless plug: http://yade.berlios.de/ ).
> > So I'm comparing your work against mine :)
> 
> Wow. It looks impressive. I note the intention in the schedule; Milestone 4 to
> "Use physical units, from dimnum or SIunits library".
> What is the status of that? Naturally having spent some time on proposed Boost
> Pqs library  http://tinyurl.com/7m5l8I would be interested to see if it might be
> used for this purpose.

well no progress on the front of physical units :) I was searching for
some useful physical units library and found those two. I also found
your pqs library, but (a year ago) after a brief look at it I decided
against it. But now that's the history :) When I have a direct contact
with the library author (you), and because Pqs will be a part of boost -
suddenly Pqs has huge advantage over other libraries, so I will try to
use it :) But not now, it is a future target ;) That URL
http://tinyurl.com/7m5l8I is not working....

 
> > 1. in file vect_def.hpp you have default constructor:
>
> How about adding a compile time switch (macro) to enable/disable zero
> initialisation?. That could also be used to measure the effect of initialisation
> on performance too.

I like that idea very much.

> > 2. in file vect.hpp, function magnitude( vect<Value_type> const & v), I
> > think you should add std:: there
> >
> > Value_type result = std::sqrt(v.x * v.x + v.y * v.y + v.z * v.z);
> 
> The intent is to allow the correct sqrt to be selected by ADL. The quantities in
> pqs have a sqrt function in boost::pqs namespace. It might be worth adding using
> std::sqrt though.

hmm, now I remember that I've seen in some boost sources STD with
capital letters, for example look at boost/serialization/vector.hpp,
there are those lines:

#if defined(__SGI_STL_PORT) || defined(_STLPORT_VERSION)
#define STD _STLP_STD
#else
#define STD std
#endif

and later STD:: is used. Maybe we should do the same...

 
> > 3. speaking of length, a function that returns squared length is very
> > useful too, because it is faster. And people often use squared length when
> > they just compare to get something longer (or shorter) but don't
> > care about the actual length.
> 
> Right. Is that the same as what norm(std::complex) does ?

yes :) I didn't know that, and from STL manual (confirmed by looking at
the source code):

     norm
     template<class Ty>
          Ty norm(const complex<Ty>& left);
     The function returns the squared magnitude of left.

this also means that name "magnitude" is more standard than "length"

> > 4. How about negation? it's operator-(T), it negates the argument.
> Yes I should add that.
> 
> > 5. how about operator-=() ; operator+=() ; operator*=() ; operator/=() ?
> Yes I must add those, assuming the argument is an arithmetic type.

yes, exactly :)

> > 6. maybe I don't understand something here, but why do you always
> > declare temporary variable "result" just to return it one line later? Is
> > it necessery for type checking? (function already has declared it's
> > return type).
> 
> Its an attempt to get the so-called Named Return Value Optimisation:
> (Heres first reference I found on it in quick search)
> http://blogs.msdn.com/slippman/archive/2004/02/03/66739.aspx

I've read this article, and it says that this optimisation is meant to
help a badly written code. Article says: "What's wrong with this
optimization is that you don't know if it will or will not be triggered"

There could be a problem with RVO, so better we test it later that there
is no burden of extra constructor call (boost is meant to work with
various compilers on different platforms) ...

> > 7. normalize() function is useful, it normalizes the vector and returns its
> > length:
> Ok. Is that useful for rotation axis for instance?

yes, eg. for constructing a new quaternion from axis/angle. Or also to
get a unit direction vector for a force (in Newton unit ;) to be applied
on a body.

 
> > 8. when working with volumes it is useful to quickly find encompassing
> > volume by performing scalar_max() on all volumes in question:
> >
> >    template <typename T>
> >    vect<BOOST_TYPEOF_TPL( T() )>
> >    scalar_max(
> >        const vect<T>& lhs,
> >        const vect<T>& rhs
> >    ){
> >        return vect<BOOST_TYPEOF_TPL( T() )>(
> >            std::max(lhs.x , rhs.x),
> >            std::max(lhs.y , rhs.y),
> >            std::max(lhs.z , rhs.z)
> >        );
> >    }
> 
> Ok. Is that sometimes implemented as a logical operator e.g Or is max and And is
> min?
> Mind you that may not be comprehensible for  vect...

Hm, I've never heard about using operators | and & for vectors for that
purpose. But when I think about it, it's maybe reasonable. Also in my
code snippet above I used name scalar_max(), although in yade the
function is just called max() - I just wasn't sure if that short name
would be clear.... Nevertheless this function is useful, and doesn't
matter too much how to call it, just ensure it's there. To be honest my
vote would go for a name "max()" :) Operators | and & are reasonable too
- although they may be not too quick to find when someone is browsing a
quick-reference to find what he needs. So whichever you prefer.

 
> > 9. scalar_min() can be sometimes useful too, and to be complete both
> > functions should be added, I think...
> 
> Ok.
> 
> > 10. what is that semicolon at the end of line 107 in vect.hpp? ;)), that '-'
> > in
> > line 106 looks suspicious too ;)
> 
> Yes they are typos. There need to be some tests!
> 
> > 11. unit_cross_product() is useful in openGL, it's just like
> > cross_product, but the result is normalized:
> > 12. we have a dot product and cross product, but how about diagonal
> > multiplication of two vectors?
> Ok. All those should be added..

nice :)
 

> > looks like that's all for now :) In your previous post you mentioned
> > quaternions, but I don't see quaternions in this geometryAL_01.zip file.
> > I'm not using matrices to rotate in my code, but quaternions, so (for
> > now) I don't have suggestions about rc_matrix. I'm only curious why the
> > "rc_" prefix to the name? I'm looking forward to see quaternions.hpp :))
> 
> the rc_ prefix stands for  Row-Column format (as opposed to Column-Row format).
> It might be best to use the OpenGL format  but I cant remember which format they
> use now. Graphics books Other Libs often use  the other format to OpenGL.  In
> fact maybe it would be best to have both and with corresponding row and column
> vectors.

hm, I don't use matrices too much so I don't know. But speaking about
OpenGL compatibility rings a bell in my mind about vectors. Many opengl
functions take as an argument a pointer to array, like float* or double*
or int*. For example this function: ( http://tinyurl.com/kfxwl )

void glVertex3fv( const GLfloat *v )

look at this example code:

    boost::geometry::three_d::vect vec(0.1,0.2,0.3);
    glVertex3fv(vec);              // shorter to write, faster to call
    glVertex3f(vec.x,vec.y,vec.z); // longer to write, slower to call, more code cluttering

To have this working, we need an operator T*() which returns pointer
to an array, which contains actual vector data, like this:

    template <class T>
    vect<T>::operator T* ()
    {
        return v; // [1]
    }
    template <class T>
    vect<T>::operator const T* () const
    {
        return v; // [1]
    }

[1] problem however is that we need that data "v" to be returned. Simplest
solution is to store this data internally in three_d::vect as an array:

   template <typename T>
   struct vect
   {
      /* ...*/
      T v[3];
   }

or, to have best of both worlds:

   template <typename T>
   struct vect
   {
      /* ...*/
       union  // replace line 35 in vect_def.hpp with those 5 lines, and everything will still commpile
       {
          struct { T x,y,z; };
          T v[3];
       };
   }

this anonymous union is an idea I took from one of books by Meyers or
Sutter (I cannot find it now :/ ). But Rodolfo Lima rightfully noted that:

>> The problem with this approach is that Type cannot have constructors
>> or destructors to get into an union. That's a little too restrictive, I
>> think.

and if we want to make a generally usable small vector we shouldn't have
this restriction - because someone may decide to use a custom class to
store numbers (for example to acheive greater precision). Here template
specializations can help us.

We can use anonymous union solution for template specializations like
float,double,int,short and leave general template that doesn't use
union, but just x y z. With this design operator T*() would be different for
general template (either can cast T to double, build an array and return
a pointer to it (requires extra fields in the struct vect {}, to store
converted values), or can simply return a compiler time error). And the
specialized versions of operator T*() would return a pointer to already
existing data "v" that contains fields x y z.

Other lazy solution is just to always build that array to which pointer
is returned (but then we lose whole speed advantage).

To have matrices similarly compatibile they also should be stored like
that in memory, look at function  ( http://tinyurl.com/gnfq8 ):

void glLoadMatrixd( const GLdouble *m );

as noted on that webpage:
   m - A pointer to a 4x4 matrix stored in column-major order as 16 consecutive values.

I see that boost::geometry::three_d::rc_matrix_base uses boost::tuple to
store data. To have rc_matrix compatilibe with glLoadMatrixd(), we also
need a similar operator double*(), which returns pointer to that matrix.
Then we can use template specializations for this operator. A general
template would generate a compile time error, while specializations for
float,int,double,etc could return pointer an array.

But that would require big changes in the rc_matrix code, so perhaps it
should be reconsidered whether those matrices should be
opengl-compatibile or not. (I've never used matrix for opengl, so I
won't use that...)

But certainly I'd need vector to be opengl-compatibile.

 
> As far as quaternions go I havent used them. I have used matrices firstly
> because a matrix can encompass a complicated transform such as rotation of an
> object around an arbitrary line not passing through the origin, whereas AFAIK a
> quaternion is limited to rotation only, so the object must first be translated
> so some part of the line is at the origin, then rotated, then translated back,
> whereas a matrix can express this transform directly in one single object. Is
> this the case or have I got it wrong?  however I am aware that quaternions are
> meant to be much quicker.

you are right, quaternions are fast when it comes to rotations, but all
rotations go around the origin. So sometimes translation is necessary.
Also if you haven't used them, it means that I'll have to push for
quaterions to get them in. This geometry library without quaternions
will be pretty useless for me :) All the fuss is about multiplication of
vector by quaternion.

But still I don't know how to put them in :) If there is a namespace
two_d and three_d, then what? Put them in three_d or four_d (???).
Should I extend existing boost::quaternion class by inheriting from it,
and adding necessary function calls, or just modify it, hoping that the
author of boost::quaternion will accept my modifications?

> > I'll try to switch my program (yade is over 30000 lines of C++ ;) into
> > this library, as a part of helping in getting this library done.
> 
> Wow that sounds like a lot of work. It sounds like your own graphics library is
> sophisticated,  so it will be interesting to see some of that and
> perhaps the work should consist more of incorporating others efforts into yours 
> ;-)

oh, it's not graphics :) Actually opengl part is not too sophisticated,
just displaying some opengl primitives. Most important part is where
physical calculations are done - collisions of bodies, calculating
forces, geometry overlapping, etc... That's where those three_d::vectors
and quaterions are used.

sorry that this post got so incredibly long ;)
 
-- 
Janek Kozicki                                                         |