[Dune] Intersection area and corners
Marco Cisternino
marco.cisternino at optimad.it
Wed Jun 1 16:01:21 CEST 2016
Hi Martin,
I have a restart procedure, but I try to give the simplest version of the code reproducing the error.
However, I'll try to rebuild dune with assertion in order to see if even with dune 2.3 I have the same assertion interruption.
I'll try to move to 2.4, but I need to say if the dune-subgrid di compatible with this version, because I need it.
Did you need to change a lot of things in my code to run with dune 2.4?
What do you mean with bugs in initialization order??
Thanks a lot for your help.
Marco
________________________________
Da: Martin Nolte <nolte at mathematik.uni-freiburg.de>
Inviato: mercoledì 1 giugno 2016 12.29
A: Dune
Oggetto: Re: [Dune] Intersection area and corners
Hi Macro,
I ported you code to the master branch (which was not too painful). For me,
the program results in an assertion (see output below), so I cannot reproduce
your bug.
Unfortunately, your "small example" is really huge (> 4000 lines of code).
This is quite a problem, because I have to understand (and trust) your example
code, before I can start hunting bugs in GeometryGrid. At the moment, I do
neither understand the code, nor do I trust it (my compiler spat out a ton of
warnings, including bugs in initialization order). It also seems to me that
the code can be further reduced.
Another problem is that the code runs quite some time in debug mode (to
compute some parametrizations). Couldn't they be precompiled offline (on your
machine) and be shipped as data?
To be honest, I am a bit at a loss how to help you. We're not supporting 2.3
anymore (since 2.4 is released) and I cannot reproduce your bug in the master
branch. Instead I get a curious assertion, which indicates you are
dereferencing an end intersection iterator.
Best,
Martin
PS: This is the output I get:
Hello World! This is intersectionBug.
I am rank 0 of 1 processes!
Starting folder selection...
The folder exists. With write permissions.
The folder exists. With write permissions.
End of folder selection...
The number of element parameters is: 0
The number of vertex parameters is: 0
Has boundary parameters? 1
4018 vertices read.
0 is the size of vtxParams_
1920 elements read: 0 simplices, 0 pyramids, 0 prisms, 1920 cubes.
0 is the size of elParams_
WARNING (ignored): Could not open file 'alugrid.cfg', using default values 0 <
[balance] < 1.2, partitioning method 'ALUGRID_SpaceFillingCurve(9)'.
You are using DUNE-ALUGrid, please don't forget to cite the paper:
Alkaemper, Dedner, Kloefkorn, Nolte. The DUNE-ALUGrid Module, 2016.
Created parallel ALUGrid<3,3,cube,nonconforming> from macro grid file
'testgrid2'.
Computing parametrizations...
End of parametrizations computation...
intersectionBug:
/home/nolte/numerics/master/dune-alugrid/dune/alugrid/3d/topology.hh:201:
static int Dune::ElementTopologyMapping<type>::dune2aluFace(int) [with
Dune::ALU3dGridElementType type = (Dune::ALU3dGridElementType)7u]: Assertion
`index >= 0 && index < numFaces' failed.
Aborted
On 05/30/2016 07:14 PM, Marco Cisternino wrote:
> Hi Martin,
>
> I found where the error is produced:
>
>
> it is exaclty in dune/grid/geometrygrid/cornerstorage.hh in the
> IntersectionCoordVector method
>
> template< std::size_t size >
> void calculate ( array< Coordinate, size > (&corners) ) const
>
> Here corners are computed by using
>
> corners[ i ] = elementGeometry_.global( hostLocalGeometry_.corner( i ) )
>
>
> hostLocalGeometry_.corner( i ) gives the right coordinates in reference
> element of corner i of the intersection, but the global() method gives the
> wrong result
>
> In global() method this statement is used to compute the return
> mapping_->global( local )
>
> if the mapping is affine (and this is the case) the method starts from
> physical position of 0-corner of the element and compute the position of the
> current corner by
> jacobianTransposed_.umtv( local, global )
>
> this is where the error is made: local is right and input global too, but
> for only one element, namely element 43 and for only 2 corner, namely
> {0,1,1}, i.e. corner 6 ,and {1,1,1}, i.e. corner 7 the results is wrong
> I would say wrong jacobian but it works fine for all the other 6 corner of the
> element (and the algorithm to compute it works for all the other intersections
> in the grid!!)
>
> Let me remark that elementGeometry_ member in IntersectionCoordVector has the
> right physical corners relatives to {0,1,1}, i.e. corner 6 ,and {1,1,1}, i.e.
> corner 7. So, it is a little bit weird to recompute intersection corners
> starting from corner 0, instead of using corners stored in elementGeometry_.
> Could you explain me why you chose this way?
>
> Thanks for any help!
> Bests,
>
> Marco
>
>
> ------------------------------------------------------------------------------
> *Da:* Marco Cisternino <marco.cisternino at optimad.it>
> *Inviato:* venerdì 27 maggio 2016 19.13
> *A:* Martin Nolte; Dune
> *Oggetto:* Re: [Dune] Intersection area and corners
>
>
> Hi Martin,
>
> I prepared the smallest case I could do.
>
> It is a bit slow because of the computation of the element parametrizations.
>
> It will print at stdout the corners of the element (index 43) having the sick
> intersection from the geometry grid and from the host one and the corners of
> the sick intersection from both grids.
>
> The intersection has index 3 in this element.
>
> As you can see in the geometry grid case corners 0 and 1 of the intersection
> are not at the same position of corners 6 and 7 of the element.
>
> In the case of the host grid everything works fine.
>
> I would no say the bug is in dune, but I cannot understand why the
> intersection is not coherent with its element.
>
> It is not evident from this code but the neighbour through this intersection
> is not coherent too. Its index is 1025.
>
>
> I remark that using the leaf mapper of the geometry grid crashes the run
> giving segmantation fault, see line 142 of intersectionBug.cc.
> All the test I made are serial.
>
> If you need more feel free to ask.
> Thanks for your help,
>
> Marco
> ------------------------------------------------------------------------------
> *Da:* Marco Cisternino <marco.cisternino at optimad.it>
> *Inviato:* giovedì 26 maggio 2016 14.12
> *A:* Martin Nolte; Dune
> *Oggetto:* Re: [Dune] Intersection area and corners
>
>
> Hi Martin,
>
> I'll try to set up a small example.
>
> The weird thing is that I have 2k elements in my coarse grid and this problem
> occurs for only 1 intersection. That's why I would like to see how the
> intersection is built, in order to understand if my code wrongly changes that
> intersection object.
>
> Thanks.
>
>
> Marco
>
>
>
> ------------------------------------------------------------------------------
> *Da:* Martin Nolte <nolte at mathematik.uni-freiburg.de>
> *Inviato:* mercoledì 25 maggio 2016 17.47
> *A:* Dune
> *Oggetto:* Re: [Dune] Intersection area and corners
>
> Hi Marco,
>
> the intersection of GeometryGrid is based on the intersection of the host
> grid. Then, you use the geometryInInside and the element geometry (of the
> inside element) to build the global geometry of the intersection. So, in a
> conforming situation, the corners of the intersection should be a subset of
> the corners of the inside element.
>
> Therefore, I consider your question a bug. Can you provide a small example?
>
> Best,
>
> Martin
>
> On 05/25/2016 03:59 PM, Marco Cisternino wrote:
>> Thank you Martin for your quick reply.
>>
>> No, for the configuration I'm debugging I'm not dealing with non-conforming
>> elements. That's why I cannot really understand why intersection and element
>> are not coherent.
>>
>> How the intersection is built? I could look in the constructor to see why its
>> corners are not among the element corners.
>>
>> Thanks again
>>
>>
>> Marco
>>
>>
>>
>>
>> ------------------------------------------------------------------------------
>> *Da:* Martin Nolte <nolte at mathematik.uni-freiburg.de>
>> *Inviato:* mercoledì 25 maggio 2016 11.00
>> *A:* Dune
>> *Oggetto:* Re: [Dune] Intersection area and corners
>>
>> Hi Marco,
>>
>> I am not sure whether I get your question correctly. You seem to be comparing
>> the geometry of the inside element to the geometry of the intersection.
>>
>> This is only reasonable, if the intersection is conforming. Otherwise I expect
>> the corners to be different. Is it possible, that your intersection is
>> non-conforming?
>>
>> Please note the using GeometryGrid with non-conforming grid views is a bit
>> problematic. Basically you are creating a first order Lagrange function over a
>> non-conforming grid. If you don't constrain hanging nodes, the position of a
>> hanging node might not be on the edge of the bigger element but slightly moved.
>>
>> You might also want to have a look at the following merge request:
>>
>> https://gitlab.dune-project.org/core/dune-grid/merge_requests/60
[https://gitlab.dune-project.org/assets/gitlab_logo-cdf021b35c4e6bb149e26460f26fae81e80552bc879179dd80e9e9266b14e894.png]<https://gitlab.dune-project.org/core/dune-grid/merge_requests/60>
Bugfix/geometrygrid various (!60) · Merge Requests · Core Modules / dune-grid<https://gitlab.dune-project.org/core/dune-grid/merge_requests/60>
gitlab.dune-project.org
Grid Interface and Implementations
> <https://gitlab.dune-project.org/core/dune-grid/merge_requests/60>
>
> Bugfix/geometrygrid various (!60) · Merge Requests · Core Modules / dune-grid
> <https://gitlab.dune-project.org/core/dune-grid/merge_requests/60>
> gitlab.dune-project.org
> Grid Interface and Implementations
>
>
>
>> <https://gitlab.dune-project.org/core/dune-grid/merge_requests/60>
>>
>> Bugfix/geometrygrid various (!60) · Merge Requests · Core Modules / dune-grid
>> <https://gitlab.dune-project.org/core/dune-grid/merge_requests/60>
>> gitlab.dune-project.org
>> Grid Interface and Implementations
>>
>>
>>
>>
>> It deals with normals in non-conforming grids, which might also be buggy in
>> DUNE 2.3.
>>
>> Best,
>>
>> Martin
>>
>> On 05/25/2016 10:44 AM, Marco Cisternino wrote:
>>> Good morning,
>>> I have a problem with one intersection in my dune application.
>>> The area I get is not what I expect.
>>> So I tried to look inside the Intersection object using gdb
>>> It is important to say that my grid is a GeometryGrid with a DiscreteFunction to parametrize the elements.
>>>
>>> If I look at these two things:
>>>
>>> - *(in.real.insideGeo_.mapping_)
>>> - *(iGeo.realGeometry.mapping_)
>>>
>>> where the in variable type is
>>> GeometryGrid::LeafGridView::IntersectionIterator::Intersection
>>> and iGeo is
>>> GeometryGrid::LeafGridView::IntersectionIterator::Intersection::Geometry
>>>
>>> I see different corners stored in. I mean two of the 4 corners given by *(iGeo.realGeometry.mapping_) are different from the two correspondant corners given by *(in.real.insideGeo_.mapping_). All the other corners are the same.
>>>
>>> If I compute the intersection area using corners given by *(in.real.insideGeo_.mapping_), I get the right area. On the other hand, if I compute the intersection area using the corners given by *(iGeo.realGeometry.mapping_), the area is wrong.
>>>
>>> I want to specify that the corners given by *(in.real.insideGeo_.mapping_) are the same I get using e.geometry().corners(i) where e is of type
>>> GeometryGrid::LeafGridView::template Codim<0>::template Partition<Dune::InteriorBorder_Partition>::Iterator::Entity
>>>
>>> Could anyone help me please? Am I doing something wrong?
>>> Ask for more details if you need.
>>>
>>> I'm using DUNE 2.3 with ALUGrid 1.52. I know it is old, but at the moment I cannot migrate to the new Dune releases.
>>>
>>> Thanks a lot!
>>>
>>> Marco
>>>
>>> _______________________________________________
>>> Dune mailing list
>>> Dune at dune-project.org
>>> http://lists.dune-project.org/mailman/listinfo/dune
>>>
>>
>> --
>> Dr. Martin Nolte <nolte at mathematik.uni-freiburg.de>
>>
>> Universität Freiburg phone: +49-761-203-5630
>> Abteilung für angewandte Mathematik fax: +49-761-203-5632
>> Hermann-Herder-Straße 10
>> 79104 Freiburg, Germany
>>
>>
>
> --
> Dr. Martin Nolte <nolte at mathematik.uni-freiburg.de>
>
> Universität Freiburg phone: +49-761-203-5630
> Abteilung für angewandte Mathematik fax: +49-761-203-5632
> Hermann-Herder-Straße 10
> 79104 Freiburg, Germany
>
>
--
Dr. Martin Nolte <nolte at mathematik.uni-freiburg.de>
Universität Freiburg phone: +49-761-203-5630
Abteilung für angewandte Mathematik fax: +49-761-203-5632
Hermann-Herder-Straße 10
79104 Freiburg, Germany
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <https://lists.dune-project.org/pipermail/dune/attachments/20160601/5f58e8e2/attachment.htm>
More information about the Dune
mailing list