Author Topic: Add tolerance as a default argument for IsInside method  (Read 568 times)

Ignasi de Pouplana

  • Newbie
  • *
  • Posts: 21
Add tolerance as a default argument for IsInside method
« on: January 02, 2017, 04:34:23 PM »
Hello everybody,

I was using the method IsInside and I was wondering whether it would be a good idea to replace the current hard-coded tolerances by a default argument. For example, in the case of the triangle_2d_3.h now we have:

Code: [Select]
    virtual bool IsInside( const CoordinatesArrayType& rPoint, CoordinatesArrayType& rResult )
    {
        const double zero = 1E-14;
        this->PointLocalCoordinates( rResult, rPoint );
        if( ( rResult[0] >= (0.0-zero) ) && ( rResult[0] <= 1.0 + zero ) )
            if( ( rResult[1] >= 0.0-zero ) && (rResult[1] <= 1.0 + zero ) )
                if(((1.0-(rResult[0] + rResult[1])) >= 0.0-zero) &&  ((1.0-(rResult[0] + rResult[1])) <= 1.0 + zero))
                    return true;

        return false;
    }

and what I propose is:

Code: [Select]
    virtual bool IsInside( const CoordinatesArrayType& rPoint, CoordinatesArrayType& rResult , const double zero=1.0e-14)
    {
        this->PointLocalCoordinates( rResult, rPoint );
        if( ( rResult[0] >= (0.0-zero) ) && ( rResult[0] <= 1.0 + zero ) )
            if( ( rResult[1] >= 0.0-zero ) && (rResult[1] <= 1.0 + zero ) )
                if(((1.0-(rResult[0] + rResult[1])) >= 0.0-zero) &&  ((1.0-(rResult[0] + rResult[1])) <= 1.0 + zero))
                    return true;

        return false;
    }

I have not committed anything yet, but if you think there is no inconvenience with this change, I would like to commit it (for every geometry). 

I'll be waiting to hear your opinions.

Regards and happy new year,

Ignasi

riccardo

  • Global Moderator
  • Newbie
  • *****
  • Posts: 47
Re: Add tolerance as a default argument for IsInside method
« Reply #1 on: January 05, 2017, 05:11:19 PM »
hi Ignaci,

to me, your proposal makes sense, so you have my +1

however committing it would imply taking care of changing the geometry base class and changing at once all of the geometries. Without this it would break the current implementation.

Pooyan, Jordi, your opinion?

cheers
Riccardo

Ignasi de Pouplana

  • Newbie
  • *
  • Posts: 21
Re: Add tolerance as a default argument for IsInside method
« Reply #2 on: January 10, 2017, 10:36:24 AM »
Ok Riccardo,

if Pooyan and Jordi find it appropriate I will change all the geometries along with the base class (I have it already prepared).

Cheers,

Ignasi

pooyan

  • Global Moderator
  • Newbie
  • *****
  • Posts: 33
Re: Add tolerance as a default argument for IsInside method
« Reply #3 on: January 12, 2017, 11:38:49 AM »
I agree with the idea. for the implementation I suggest to use std::numeric_limits<double>::epsilon() instead of 1e-14 as default value and name the variable as Tolerance (not zero!)

Bests,

Ignasi de Pouplana

  • Newbie
  • *
  • Posts: 21
Re: Add tolerance as a default argument for IsInside method
« Reply #4 on: January 12, 2017, 12:38:24 PM »
I'll take into account your suggestions and after checking that everything compiles (both in linux and windows) I'll commit the changes.

Regards.

jcotela

  • Newbie
  • *
  • Posts: 4
Re: Add tolerance as a default argument for IsInside method
« Reply #5 on: January 30, 2017, 05:50:17 PM »
Hi Ignasi,
while I agree with the proposal in principle, I'm finding out that what is currently on svn is not working as intended. If I'm interpreting my results correctly, the problem is as follows:

- For anything except linear triangles/tetrahedra, determining the local coordinates is an iterative process: the PointLocalCoordinates function performs an iterative solution (currently implemented with a hard-coded tolerance of 1e-8).
- Then, when I use the IsInside function for tetrahedra without arguments, "good" PointLocalCoordinates results (to a tolerance of 1e-8) fail IsInside (to machine tolerance) and, as a result, the routine rejects the point as being outside.

I'm not sure how we want to solve this. An obvious solution would be to pass the IsInside tolerance to PointLocalCoordinates, but:
- Passing a tolerance to PointLocalCoordinates makes no sense for elements where it is exact (triangles, tetrahedra).
- I'm a bit reluctant to solve the Local Coordinates problem to machine tolerance by default (I fear it will be unnecessarily expensive, or that convergence might be prevented by cumulative errors for complex elements).

Does anyone else have an opinion on this?

Ignasi de Pouplana

  • Newbie
  • *
  • Posts: 21
Re: Add tolerance as a default argument for IsInside method
« Reply #6 on: January 30, 2017, 07:04:49 PM »
I agree with Jordi.

Using std::numeric_limits<double>::epsilon() as the default tolerance can be incompatible with the hard-coded tolerance of the PointLocalCoordinates method...

As a means to prevent problems as the one reported by Jordi, I would set 1.0e-08 as the default tolerance for all the geometries and, if someone wants a smaller tolerance, just call IsInside with it in the argument.

What do you think?