Dear Giovanni, Thank you very much for the clear and simple test case! I've played around with it and I don't see any fault in the way that you're trying to use the mapping functionality, so I think that its quite clear that there's something wrong with Mapping::transform_real_to_unit_cell. I've tested this with both deal.II 8.5 and the development version and the problem exists in both cases. I've opened a issue on GitHub <https://github.com/dealii/dealii/issues/4654> that uses your example to demonstrate the problem. I'll report back here when we make some headway towards solving the issue!
Best regards, Jean-Paul On Tuesday, July 25, 2017 at 12:29:51 AM UTC+2, Giovanni Di Ilio wrote: > > Sorry, I thought I worked around the problem but actually I did not. The > mapping fails again when I work with a refined mesh, for both the cases of > a mesh read in from gmsh and a mesh generated with Deal. > > I try to explain better what I mean: > 1) I generate a grid, using: > *GridGenerator::subdivided_hyper_rectangle(triangulation,repetitions,p1,p2)* > ; > 2) I refine one time the grid with: *triangulation.refine_global(1)* ; > 3) I run over the cells and I try to get the mapped position of a given > point, which I set to be in the centre of the considered cell: > *mapping.transform_real_to_unit_cell*. > > Since the point I am testing is exactly in the centre of the cell I expect > the mapped position, in the unit cell, to be x = 0.5; y = 0.5. > However, for some "random" cell, again, the mapping fails. > > Here I am not dealing anymore with the digits of accuracy. > The numbers I am dealing with are the coordinates of the vertexes provided > by the GridGenerator and the coordinates of a testing point, which have > just a couple of digits after comma. Moreover, the spacing between vertexes > for the cells I am considering is of the order of 0.1. > > I attach a very simple example, where I consider a uniform grid of 5 x 5 > elements (repetitions[0]=repetitions[1]=5), refined one time, so 100 > elements in total. > For cell with index 96 of my example, which vertexes are: > vertex(0): 0.8 0.8 > vertex(1): 0.9 0.8 > vertex(2): 0.8 0.9 > vertex(3): 0.9 0.9 > I test the point (0.85, 0.85) and the mapped position in the unit cell > results to be: (0.5, 0.46875) instead of (0.5, 0.5). > > You may notice that, by running the example with refineglobal(false), the > mapping on an equivalent grid made of 10 x 10 elements, with no refinements > (repetitions[0]=repetitions[1]=10), works just fine. > I don't see any difference between the two cases. > What is missing here? Can you help me please? > > Thank you very much, > Best, > > Giovanni > > > > On Saturday, July 22, 2017 at 8:57:57 PM UTC+2, Wolfgang Bangerth wrote: >> >> On 07/22/2017 04:42 AM, Giovanni Di Ilio wrote: >> > >> > Yes, I agree, the issue is definitely related to round-off error. But I >> still >> > don't understand why the mapping fails, if the vertex coordinates read >> in from >> > my .msh file are apparently stored as doubles... >> > Anyway, yes, the problem is easily overcome by "manipulating" the >> values >> > before mapping. >> >> The numbers you consider differ in the last two or three bits when stored >> as >> double precision. Any computations you do with them, for example in the >> mapping, will then necessarily be completely off. >> >> Best >> W. >> >> -- >> ------------------------------------------------------------------------ >> Wolfgang Bangerth email: bang...@colostate.edu >> www: http://www.math.colostate.edu/~bangerth/ >> >> -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.