Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[Jean-Paul Pelteret]

Seconded! It’s a really nice result. I’m glad that you managed to solve your 
problems, Bruno.

Best,
Jean-Paul

> On 09 Apr 2020, at 11:07, luca.heltai  wrote:
> 
> That’s fanstastic!
> 
> :)
> 
> You sould put this image in the image gallery!
> 
> Luca.
> 
>> On 8 Apr 2020, at 16:29, Bruno Blais  wrote:
>> 
>> Dear all,
>> Thank you very much for the help and support. I have managed to make it work 
>> in 3D within Lethe.
>> A quick image of the application of CAD mesh adaptation using multiple CADs 
>> on a complex mixer :)!
>> Thank  you for your help (especially you Luca :)!)
>> 
>> 
>> 
>> 
>> On Sunday, 5 April 2020 12:05:17 UTC-4, Bruno Blais wrote:
>> Dear Jean-Paul, what you wrote in the FAQ is very clear. If I find 
>> additional elements as I play with the functionnality I will add it to the 
>> FAQ. Thanks for adding this section.
>> Best
>> Bruno
>> 
>> 
>> On Friday, 3 April 2020 16:11:05 UTC-4, Jean-Paul Pelteret wrote:
>> Dear Bruno,
>> 
>> I guess that it would be best to mention this in the tutorial itself, but in 
>> the mean time I wrote a segment about it in the FAQ:
>> https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#how-do-i-ensure-that-my-refined-grid-is-correct-when-using-cad-geometry
>> Please feel free to edit it as you see fit. Never having used this 
>> functionality in the library before, I hope that I have correctly 
>> interpreted and massages what Luca had detailed.
>> 
>> Best,
>> Jean-Paul
>> 
>>> On 03 Apr 2020, at 15:55, Bruno Blais  wrote:
>>> 
>>> Dear Jean-Paul,
>>> Do you think it should be added to the wiki page or as additional 
>>> information to step-54? I have a bit on my hands, but i'd be glad to work 
>>> to make the information readily available. Clearly I would not have found 
>>> that out without Luca's help. 
>>> Best
>>> Bruno
>>> 
>>> 
>>> On Friday, 3 April 2020 09:43:33 UTC-4, Jean-Paul Pelteret wrote:
>>> Dear Bruno,
>>> 
>>> I’m glad that you managed to solve this issue with Luca’s aid. I think that 
>>> the information and workflow that Luca gave was really valuable. Would you 
>>> care to add this to the Wiki page?
>>> 
>>> Best,
>>> Jean-Paul
>>> 
 On 03 Apr 2020, at 06:03, Bruno Blais  wrote:
 
 Ciao Luca,
 It works perfectly now even with regular GMSH triangulation (see movie). 
 Now I understand the difference and how to set it correctly.
 Thank you very much for the help :)!
 
 
 On Thursday, 2 April 2020 11:31:42 UTC-4, luca.heltai wrote:
 Bruno, it seems like you are attaching your faces to the *boundary* 
 manifold. Let me try to explain what is happening. 
 
 In the code for step-54, there is a wire that is used to describe the 
 manifold of the *boundary* of the surface (a curve of dimension one 
 embedded in dimension three). This is generated by the wire that 
 identifies the boundary of the TopoDS shape of the face. The manifold 
 attached to it is an ArcLengthProjectionManifold, where mid points on the 
 curve are added by computing the distance in arc length, and taking the 
 point in the middle w.r.t. this length. 
 
 In your code with more than one cell, you are not assigning *just the 
 boundary faces* to the wire, but *all* faces (including the one in the 
 middle of your surface, delimitating the two cells). Intermediate points 
 on those lines will be projected (correctly) on the boundary wire (that’s 
 what your movie shows). 
 
 Of course any point that is in the middle of the surface does *not* belong 
 to the *boundary wire*, and you’ll get an exception when trying to 
 construct the midpoint of an *edge* that belongs to the interior of the 
 domain, with manifold id 2 (which corresponds to a manifold that describes 
 *the boundary* of your topology, not its surface). 
 
 With two cells, the internal edge (the only edge that does not belong to 
 the boundary) should *not* have id 2. If you set its id to 2, the result 
 is that a point which is in the middle of the edge, is actually in the 
 midlle of the two points *when running along the curve that describes the 
 boundary*, hence you get the distortion you show in the movie. 
 
 Try adding a check to the faces in which you set the manifold id to 2 only 
 if the face is at the boundary. 
 
 The principle is: 
 1. start from the lowest codimension objects, identify how to deform them. 
 In your case, cells of a Tria<2,3> are quads, that should deform as a 
 TopoDS_FACE. 
 2. Attach your favorite manifold to the TopoDS_Face (I’d personally only 
 use NormalToMeshProjectionManifold) using 
 cell->set_all_manifold_ids(manifold_id) (Notice the use of 
 set_all_manifold_ids, and **not** set_manifold_id: you want all children 
 of these objects to belong to the same manifold, in particular you want 
 all faces that are internal, i.e., that are shared between 

Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[luca.heltai]

That’s fanstastic!

:)

You sould put this image in the image gallery!

Luca.

> On 8 Apr 2020, at 16:29, Bruno Blais  wrote:
> 
> Dear all,
> Thank you very much for the help and support. I have managed to make it work 
> in 3D within Lethe.
> A quick image of the application of CAD mesh adaptation using multiple CADs 
> on a complex mixer :)!
> Thank  you for your help (especially you Luca :)!)
> 
> 
> 
> 
> On Sunday, 5 April 2020 12:05:17 UTC-4, Bruno Blais wrote:
> Dear Jean-Paul, what you wrote in the FAQ is very clear. If I find additional 
> elements as I play with the functionnality I will add it to the FAQ. Thanks 
> for adding this section.
> Best
> Bruno
> 
> 
> On Friday, 3 April 2020 16:11:05 UTC-4, Jean-Paul Pelteret wrote:
> Dear Bruno,
> 
> I guess that it would be best to mention this in the tutorial itself, but in 
> the mean time I wrote a segment about it in the FAQ:
> https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#how-do-i-ensure-that-my-refined-grid-is-correct-when-using-cad-geometry
> Please feel free to edit it as you see fit. Never having used this 
> functionality in the library before, I hope that I have correctly interpreted 
> and massages what Luca had detailed.
> 
> Best,
> Jean-Paul
> 
>> On 03 Apr 2020, at 15:55, Bruno Blais  wrote:
>> 
>> Dear Jean-Paul,
>> Do you think it should be added to the wiki page or as additional 
>> information to step-54? I have a bit on my hands, but i'd be glad to work to 
>> make the information readily available. Clearly I would not have found that 
>> out without Luca's help. 
>> Best
>> Bruno
>> 
>> 
>> On Friday, 3 April 2020 09:43:33 UTC-4, Jean-Paul Pelteret wrote:
>> Dear Bruno,
>> 
>> I’m glad that you managed to solve this issue with Luca’s aid. I think that 
>> the information and workflow that Luca gave was really valuable. Would you 
>> care to add this to the Wiki page?
>> 
>> Best,
>> Jean-Paul
>> 
>>> On 03 Apr 2020, at 06:03, Bruno Blais  wrote:
>>> 
>>> Ciao Luca,
>>> It works perfectly now even with regular GMSH triangulation (see movie). 
>>> Now I understand the difference and how to set it correctly.
>>> Thank you very much for the help :)!
>>> 
>>> 
>>> On Thursday, 2 April 2020 11:31:42 UTC-4, luca.heltai wrote:
>>> Bruno, it seems like you are attaching your faces to the *boundary* 
>>> manifold. Let me try to explain what is happening. 
>>> 
>>> In the code for step-54, there is a wire that is used to describe the 
>>> manifold of the *boundary* of the surface (a curve of dimension one 
>>> embedded in dimension three). This is generated by the wire that identifies 
>>> the boundary of the TopoDS shape of the face. The manifold attached to it 
>>> is an ArcLengthProjectionManifold, where mid points on the curve are added 
>>> by computing the distance in arc length, and taking the point in the middle 
>>> w.r.t. this length. 
>>> 
>>> In your code with more than one cell, you are not assigning *just the 
>>> boundary faces* to the wire, but *all* faces (including the one in the 
>>> middle of your surface, delimitating the two cells). Intermediate points on 
>>> those lines will be projected (correctly) on the boundary wire (that’s what 
>>> your movie shows). 
>>> 
>>> Of course any point that is in the middle of the surface does *not* belong 
>>> to the *boundary wire*, and you’ll get an exception when trying to 
>>> construct the midpoint of an *edge* that belongs to the interior of the 
>>> domain, with manifold id 2 (which corresponds to a manifold that describes 
>>> *the boundary* of your topology, not its surface). 
>>> 
>>> With two cells, the internal edge (the only edge that does not belong to 
>>> the boundary) should *not* have id 2. If you set its id to 2, the result is 
>>> that a point which is in the middle of the edge, is actually in the midlle 
>>> of the two points *when running along the curve that describes the 
>>> boundary*, hence you get the distortion you show in the movie. 
>>> 
>>> Try adding a check to the faces in which you set the manifold id to 2 only 
>>> if the face is at the boundary. 
>>> 
>>> The principle is: 
>>> 1. start from the lowest codimension objects, identify how to deform them. 
>>> In your case, cells of a Tria<2,3> are quads, that should deform as a 
>>> TopoDS_FACE. 
>>> 2. Attach your favorite manifold to the TopoDS_Face (I’d personally only 
>>> use NormalToMeshProjectionManifold) using 
>>> cell->set_all_manifold_ids(manifold_id) (Notice the use of 
>>> set_all_manifold_ids, and **not** set_manifold_id: you want all children of 
>>> these objects to belong to the same manifold, in particular you want all 
>>> faces that are internal, i.e., that are shared between two obects with the 
>>> same manifold id, to inherit the same manifold id). 
>>> 3. Go one codimension lower: in your case, curves (for 3d meshes, 
>>> surfaces). Follow the same rule as above, setting all_manifold_ids on faces 
>>> that you know should follow a known codimension one 

Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[Bruno Blais]

Dear Jean-Paul, what you wrote in the FAQ is very clear. If I find 
additional elements as I play with the functionnality I will add it to the 
FAQ. Thanks for adding this section.
Best
Bruno


On Friday, 3 April 2020 16:11:05 UTC-4, Jean-Paul Pelteret wrote:
>
> Dear Bruno,
>
> I guess that it would be best to mention this in the tutorial itself, but 
> in the mean time I wrote a segment about it in the FAQ:
>
> https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#how-do-i-ensure-that-my-refined-grid-is-correct-when-using-cad-geometry
> Please feel free to edit it as you see fit. Never having used this 
> functionality in the library before, I hope that I have correctly 
> interpreted and massages what Luca had detailed.
>
> Best,
> Jean-Paul
>
> On 03 Apr 2020, at 15:55, Bruno Blais > 
> wrote:
>
> Dear Jean-Paul,
> Do you think it should be added to the wiki page or as additional 
> information to step-54? I have a bit on my hands, but i'd be glad to work 
> to make the information readily available. Clearly I would not have found 
> that out without Luca's help. 
> Best
> Bruno
>
>
> On Friday, 3 April 2020 09:43:33 UTC-4, Jean-Paul Pelteret wrote:
>>
>> Dear Bruno,
>>
>> I’m glad that you managed to solve this issue with Luca’s aid. I think 
>> that the information and workflow that Luca gave was really valuable. Would 
>> you care to add this to the Wiki page?
>>
>> Best,
>> Jean-Paul
>>
>> On 03 Apr 2020, at 06:03, Bruno Blais  wrote:
>>
>> Ciao Luca,
>> It works perfectly now even with regular GMSH triangulation (see movie). 
>> Now I understand the difference and how to set it correctly.
>> Thank you very much for the help :)!
>>
>>
>> On Thursday, 2 April 2020 11:31:42 UTC-4, luca.heltai wrote:
>>>
>>> Bruno, it seems like you are attaching your faces to the *boundary* 
>>> manifold. Let me try to explain what is happening. 
>>>
>>> In the code for step-54, there is a wire that is used to describe the 
>>> manifold of the *boundary* of the surface (a curve of dimension one 
>>> embedded in dimension three). This is generated by the wire that identifies 
>>> the boundary of the TopoDS shape of the face. The manifold attached to it 
>>> is an ArcLengthProjectionManifold, where mid points on the curve are added 
>>> by computing the distance in arc length, and taking the point in the middle 
>>> w.r.t. this length. 
>>>
>>> In your code with more than one cell, you are not assigning *just the 
>>> boundary faces* to the wire, but *all* faces (including the one in the 
>>> middle of your surface, delimitating the two cells). Intermediate points on 
>>> those lines will be projected (correctly) on the boundary wire (that’s what 
>>> your movie shows). 
>>>
>>> Of course any point that is in the middle of the surface does *not* 
>>> belong to the *boundary wire*, and you’ll get an exception when trying to 
>>> construct the midpoint of an *edge* that belongs to the interior of the 
>>> domain, with manifold id 2 (which corresponds to a manifold that describes 
>>> *the boundary* of your topology, not its surface). 
>>>
>>> With two cells, the internal edge (the only edge that does not belong to 
>>> the boundary) should *not* have id 2. If you set its id to 2, the result is 
>>> that a point which is in the middle of the edge, is actually in the midlle 
>>> of the two points *when running along the curve that describes the 
>>> boundary*, hence you get the distortion you show in the movie. 
>>>
>>> Try adding a check to the faces in which you set the manifold id to 2 
>>> only if the face is at the boundary. 
>>>
>>> The principle is: 
>>> 1. start from the lowest codimension objects, identify how to deform 
>>> them. In your case, cells of a Tria<2,3> are quads, that should deform as a 
>>> TopoDS_FACE. 
>>> 2. Attach your favorite manifold to the TopoDS_Face (I’d personally only 
>>> use NormalToMeshProjectionManifold) using 
>>> cell->set_all_manifold_ids(manifold_id) (Notice the use of 
>>> set_all_manifold_ids, and **not** set_manifold_id: you want all children of 
>>> these objects to belong to the same manifold, in particular you want all 
>>> faces that are internal, i.e., that are shared between two obects with the 
>>> same manifold id, to inherit the same manifold id). 
>>> 3. Go one codimension lower: in your case, curves (for 3d meshes, 
>>> surfaces). Follow the same rule as above, setting all_manifold_ids on faces 
>>> that you know should follow a known codimension one geometry (a known 
>>> TopoDS_EDGE, or TopoDS_WIRE). In this case I’d only use 
>>> ArcLengthProjectionManifold objects, attached to the wires that identify 
>>> your geometry. 
>>>
>>> Doing things in this order guarantees that internal edges get the 
>>> correct manifold id, which, in your case, is not happening. 
>>>
>>> Best, 
>>> Luca. 
>>>
>>>
>>>
>>> > On 1 Apr 2020, at 17:32, Bruno Blais  wrote: 
>>> > 
>>> > So my investigation on this issue continues... 
>>> > It seems the core of my issues stem 

Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[Jean-Paul Pelteret]

Dear Bruno,

I guess that it would be best to mention this in the tutorial itself, but in 
the mean time I wrote a segment about it in the FAQ:
https://github.com/dealii/dealii/wiki/Frequently-Asked-Questions#how-do-i-ensure-that-my-refined-grid-is-correct-when-using-cad-geometry
 

Please feel free to edit it as you see fit. Never having used this 
functionality in the library before, I hope that I have correctly interpreted 
and massages what Luca had detailed.

Best,
Jean-Paul

> On 03 Apr 2020, at 15:55, Bruno Blais  wrote:
> 
> Dear Jean-Paul,
> Do you think it should be added to the wiki page or as additional information 
> to step-54? I have a bit on my hands, but i'd be glad to work to make the 
> information readily available. Clearly I would not have found that out 
> without Luca's help. 
> Best
> Bruno
> 
> 
> On Friday, 3 April 2020 09:43:33 UTC-4, Jean-Paul Pelteret wrote:
> Dear Bruno,
> 
> I’m glad that you managed to solve this issue with Luca’s aid. I think that 
> the information and workflow that Luca gave was really valuable. Would you 
> care to add this to the Wiki page?
> 
> Best,
> Jean-Paul
> 
>> On 03 Apr 2020, at 06:03, Bruno Blais > 
>> wrote:
>> 
>> Ciao Luca,
>> It works perfectly now even with regular GMSH triangulation (see movie). Now 
>> I understand the difference and how to set it correctly.
>> Thank you very much for the help :)!
>> 
>> 
>> On Thursday, 2 April 2020 11:31:42 UTC-4, luca.heltai wrote:
>> Bruno, it seems like you are attaching your faces to the *boundary* 
>> manifold. Let me try to explain what is happening. 
>> 
>> In the code for step-54, there is a wire that is used to describe the 
>> manifold of the *boundary* of the surface (a curve of dimension one embedded 
>> in dimension three). This is generated by the wire that identifies the 
>> boundary of the TopoDS shape of the face. The manifold attached to it is an 
>> ArcLengthProjectionManifold, where mid points on the curve are added by 
>> computing the distance in arc length, and taking the point in the middle 
>> w.r.t. this length. 
>> 
>> In your code with more than one cell, you are not assigning *just the 
>> boundary faces* to the wire, but *all* faces (including the one in the 
>> middle of your surface, delimitating the two cells). Intermediate points on 
>> those lines will be projected (correctly) on the boundary wire (that’s what 
>> your movie shows). 
>> 
>> Of course any point that is in the middle of the surface does *not* belong 
>> to the *boundary wire*, and you’ll get an exception when trying to construct 
>> the midpoint of an *edge* that belongs to the interior of the domain, with 
>> manifold id 2 (which corresponds to a manifold that describes *the boundary* 
>> of your topology, not its surface). 
>> 
>> With two cells, the internal edge (the only edge that does not belong to the 
>> boundary) should *not* have id 2. If you set its id to 2, the result is that 
>> a point which is in the middle of the edge, is actually in the midlle of the 
>> two points *when running along the curve that describes the boundary*, hence 
>> you get the distortion you show in the movie. 
>> 
>> Try adding a check to the faces in which you set the manifold id to 2 only 
>> if the face is at the boundary. 
>> 
>> The principle is: 
>> 1. start from the lowest codimension objects, identify how to deform them. 
>> In your case, cells of a Tria<2,3> are quads, that should deform as a 
>> TopoDS_FACE. 
>> 2. Attach your favorite manifold to the TopoDS_Face (I’d personally only use 
>> NormalToMeshProjectionManifold) using 
>> cell->set_all_manifold_ids(manifold_id) (Notice the use of 
>> set_all_manifold_ids, and **not** set_manifold_id: you want all children of 
>> these objects to belong to the same manifold, in particular you want all 
>> faces that are internal, i.e., that are shared between two obects with the 
>> same manifold id, to inherit the same manifold id). 
>> 3. Go one codimension lower: in your case, curves (for 3d meshes, surfaces). 
>> Follow the same rule as above, setting all_manifold_ids on faces that you 
>> know should follow a known codimension one geometry (a known TopoDS_EDGE, or 
>> TopoDS_WIRE). In this case I’d only use ArcLengthProjectionManifold objects, 
>> attached to the wires that identify your geometry. 
>> 
>> Doing things in this order guarantees that internal edges get the correct 
>> manifold id, which, in your case, is not happening. 
>> 
>> Best, 
>> Luca. 
>> 
>> 
>> 
>> > On 1 Apr 2020, at 17:32, Bruno Blais gmail.com 
>> > > wrote: 
>> > 
>> > So my investigation on this issue continues... 
>> > It seems the core of my issues stem from using a non-trivial mesh with 
>> > more than one cell as a starting point. 
>> > if I start with a 2 cell spacedim=2 dim=3 mesh such as: 
>> > 
>> > 
>> >  
>> > 
>> 

Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[Bruno Blais]

Dear Jean-Paul,
Do you think it should be added to the wiki page or as additional 
information to step-54? I have a bit on my hands, but i'd be glad to work 
to make the information readily available. Clearly I would not have found 
that out without Luca's help. 
Best
Bruno


On Friday, 3 April 2020 09:43:33 UTC-4, Jean-Paul Pelteret wrote:
>
> Dear Bruno,
>
> I’m glad that you managed to solve this issue with Luca’s aid. I think 
> that the information and workflow that Luca gave was really valuable. Would 
> you care to add this to the Wiki page?
>
> Best,
> Jean-Paul
>
> On 03 Apr 2020, at 06:03, Bruno Blais > 
> wrote:
>
> Ciao Luca,
> It works perfectly now even with regular GMSH triangulation (see movie). 
> Now I understand the difference and how to set it correctly.
> Thank you very much for the help :)!
>
>
> On Thursday, 2 April 2020 11:31:42 UTC-4, luca.heltai wrote:
>>
>> Bruno, it seems like you are attaching your faces to the *boundary* 
>> manifold. Let me try to explain what is happening. 
>>
>> In the code for step-54, there is a wire that is used to describe the 
>> manifold of the *boundary* of the surface (a curve of dimension one 
>> embedded in dimension three). This is generated by the wire that identifies 
>> the boundary of the TopoDS shape of the face. The manifold attached to it 
>> is an ArcLengthProjectionManifold, where mid points on the curve are added 
>> by computing the distance in arc length, and taking the point in the middle 
>> w.r.t. this length. 
>>
>> In your code with more than one cell, you are not assigning *just the 
>> boundary faces* to the wire, but *all* faces (including the one in the 
>> middle of your surface, delimitating the two cells). Intermediate points on 
>> those lines will be projected (correctly) on the boundary wire (that’s what 
>> your movie shows). 
>>
>> Of course any point that is in the middle of the surface does *not* 
>> belong to the *boundary wire*, and you’ll get an exception when trying to 
>> construct the midpoint of an *edge* that belongs to the interior of the 
>> domain, with manifold id 2 (which corresponds to a manifold that describes 
>> *the boundary* of your topology, not its surface). 
>>
>> With two cells, the internal edge (the only edge that does not belong to 
>> the boundary) should *not* have id 2. If you set its id to 2, the result is 
>> that a point which is in the middle of the edge, is actually in the midlle 
>> of the two points *when running along the curve that describes the 
>> boundary*, hence you get the distortion you show in the movie. 
>>
>> Try adding a check to the faces in which you set the manifold id to 2 
>> only if the face is at the boundary. 
>>
>> The principle is: 
>> 1. start from the lowest codimension objects, identify how to deform 
>> them. In your case, cells of a Tria<2,3> are quads, that should deform as a 
>> TopoDS_FACE. 
>> 2. Attach your favorite manifold to the TopoDS_Face (I’d personally only 
>> use NormalToMeshProjectionManifold) using 
>> cell->set_all_manifold_ids(manifold_id) (Notice the use of 
>> set_all_manifold_ids, and **not** set_manifold_id: you want all children of 
>> these objects to belong to the same manifold, in particular you want all 
>> faces that are internal, i.e., that are shared between two obects with the 
>> same manifold id, to inherit the same manifold id). 
>> 3. Go one codimension lower: in your case, curves (for 3d meshes, 
>> surfaces). Follow the same rule as above, setting all_manifold_ids on faces 
>> that you know should follow a known codimension one geometry (a known 
>> TopoDS_EDGE, or TopoDS_WIRE). In this case I’d only use 
>> ArcLengthProjectionManifold objects, attached to the wires that identify 
>> your geometry. 
>>
>> Doing things in this order guarantees that internal edges get the correct 
>> manifold id, which, in your case, is not happening. 
>>
>> Best, 
>> Luca. 
>>
>>
>>
>> > On 1 Apr 2020, at 17:32, Bruno Blais  wrote: 
>> > 
>> > So my investigation on this issue continues... 
>> > It seems the core of my issues stem from using a non-trivial mesh with 
>> more than one cell as a starting point. 
>> > if I start with a 2 cell spacedim=2 dim=3 mesh such as: 
>> > 
>> > 
>> >  
>> > 
>> > 
>> > 
>> > Which is identical to the second step of the adaptation. I get an 
>> aberrant result (see attached out_wrong.mp4). 
>> > 
>> > 
>> > 
>> > I also get issues when I try to apply the manifold on a starting msh 
>> made with GMSH that has a similar topology such as in the following image:
>>  
>> > 
>> >  
>> > 
>> > 
>> > 
>> > In this case, I get an error about a point not being on the manifold 
>> (when it clearly is)... 
>> > 
>> > So I must admit I am at loss here. Right now the code I use to set the 
>> manifold is : 
>> > 
>> > for (const auto  : tria.active_cell_iterators()) 
>> >   { 
>> > 
>> > cell->set_manifold_id(1); 
>> > 
>> > for (const auto  : cell->face_iterators()) 
>> > 
>> >   { 

Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[Jean-Paul Pelteret]

Dear Bruno,

I’m glad that you managed to solve this issue with Luca’s aid. I think that the 
information and workflow that Luca gave was really valuable. Would you care to 
add this to the Wiki page?

Best,
Jean-Paul

> On 03 Apr 2020, at 06:03, Bruno Blais  wrote:
> 
> Ciao Luca,
> It works perfectly now even with regular GMSH triangulation (see movie). Now 
> I understand the difference and how to set it correctly.
> Thank you very much for the help :)!
> 
> 
> On Thursday, 2 April 2020 11:31:42 UTC-4, luca.heltai wrote:
> Bruno, it seems like you are attaching your faces to the *boundary* manifold. 
> Let me try to explain what is happening. 
> 
> In the code for step-54, there is a wire that is used to describe the 
> manifold of the *boundary* of the surface (a curve of dimension one embedded 
> in dimension three). This is generated by the wire that identifies the 
> boundary of the TopoDS shape of the face. The manifold attached to it is an 
> ArcLengthProjectionManifold, where mid points on the curve are added by 
> computing the distance in arc length, and taking the point in the middle 
> w.r.t. this length. 
> 
> In your code with more than one cell, you are not assigning *just the 
> boundary faces* to the wire, but *all* faces (including the one in the middle 
> of your surface, delimitating the two cells). Intermediate points on those 
> lines will be projected (correctly) on the boundary wire (that’s what your 
> movie shows). 
> 
> Of course any point that is in the middle of the surface does *not* belong to 
> the *boundary wire*, and you’ll get an exception when trying to construct the 
> midpoint of an *edge* that belongs to the interior of the domain, with 
> manifold id 2 (which corresponds to a manifold that describes *the boundary* 
> of your topology, not its surface). 
> 
> With two cells, the internal edge (the only edge that does not belong to the 
> boundary) should *not* have id 2. If you set its id to 2, the result is that 
> a point which is in the middle of the edge, is actually in the midlle of the 
> two points *when running along the curve that describes the boundary*, hence 
> you get the distortion you show in the movie. 
> 
> Try adding a check to the faces in which you set the manifold id to 2 only if 
> the face is at the boundary. 
> 
> The principle is: 
> 1. start from the lowest codimension objects, identify how to deform them. In 
> your case, cells of a Tria<2,3> are quads, that should deform as a 
> TopoDS_FACE. 
> 2. Attach your favorite manifold to the TopoDS_Face (I’d personally only use 
> NormalToMeshProjectionManifold) using cell->set_all_manifold_ids(manifold_id) 
> (Notice the use of set_all_manifold_ids, and **not** set_manifold_id: you 
> want all children of these objects to belong to the same manifold, in 
> particular you want all faces that are internal, i.e., that are shared 
> between two obects with the same manifold id, to inherit the same manifold 
> id). 
> 3. Go one codimension lower: in your case, curves (for 3d meshes, surfaces). 
> Follow the same rule as above, setting all_manifold_ids on faces that you 
> know should follow a known codimension one geometry (a known TopoDS_EDGE, or 
> TopoDS_WIRE). In this case I’d only use ArcLengthProjectionManifold objects, 
> attached to the wires that identify your geometry. 
> 
> Doing things in this order guarantees that internal edges get the correct 
> manifold id, which, in your case, is not happening. 
> 
> Best, 
> Luca. 
> 
> 
> 
> > On 1 Apr 2020, at 17:32, Bruno Blais gmail.com 
> > > wrote: 
> > 
> > So my investigation on this issue continues... 
> > It seems the core of my issues stem from using a non-trivial mesh with more 
> > than one cell as a starting point. 
> > if I start with a 2 cell spacedim=2 dim=3 mesh such as: 
> > 
> > 
> >  
> > 
> > 
> > 
> > Which is identical to the second step of the adaptation. I get an aberrant 
> > result (see attached out_wrong.mp4). 
> > 
> > 
> > 
> > I also get issues when I try to apply the manifold on a starting msh made 
> > with GMSH that has a similar topology such as in the following image: 
> > 
> >  
> > 
> > 
> > 
> > In this case, I get an error about a point not being on the manifold (when 
> > it clearly is)... 
> > 
> > So I must admit I am at loss here. Right now the code I use to set the 
> > manifold is : 
> > 
> > for (const auto  : tria.active_cell_iterators()) 
> >   { 
> > 
> > cell->set_manifold_id(1); 
> > 
> > for (const auto  : cell->face_iterators()) 
> > 
> >   { 
> > 
> >   face->set_manifold_id(2); 
> > 
> >   } 
> > 
> >   } 
> > 
> > If anybody has any experience with using the OpenCascade manifold using 
> > initial meshes that have more than 1 cell, I would greatly appreciate 
> > advices :)!. 
> > 
> > Best 
> > 
> > Bruno 
> > 
> > 
> > 
> > 
> > On Wednesday, April 1, 2020 at 9:38:27 AM UTC-4, Bruno Blais wrote: 
> > 
> > So it seems part 

Re: [deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[luca.heltai]

Bruno, it seems like you are attaching your faces to the *boundary* manifold. 
Let me try to explain what is happening. 

In the code for step-54, there is a wire that is used to describe the manifold 
of the *boundary* of the surface (a curve of dimension one embedded in 
dimension three). This is generated by the wire that identifies the boundary of 
the TopoDS shape of the face. The manifold attached to it is an 
ArcLengthProjectionManifold, where mid points on the curve are added by 
computing the distance in arc length, and taking the point in the middle w.r.t. 
this length.

In your code with more than one cell, you are not assigning *just the boundary 
faces* to the wire, but *all* faces (including the one in the middle of your 
surface, delimitating the two cells). Intermediate points on those lines will 
be projected (correctly) on the boundary wire (that’s what your movie shows). 

Of course any point that is in the middle of the surface does *not* belong to 
the *boundary wire*, and you’ll get an exception when trying to construct the 
midpoint of an *edge* that belongs to the interior of the domain, with manifold 
id 2 (which corresponds to a manifold that describes *the boundary* of your 
topology, not its surface).

With two cells, the internal edge (the only edge that does not belong to the 
boundary) should *not* have id 2. If you set its id to 2, the result is that a 
point which is in the middle of the edge, is actually in the midlle of the two 
points *when running along the curve that describes the boundary*, hence you 
get the distortion you show in the movie. 

Try adding a check to the faces in which you set the manifold id to 2 only if 
the face is at the boundary.

The principle is: 
1. start from the lowest codimension objects, identify how to deform them. In 
your case, cells of a Tria<2,3> are quads, that should deform as a TopoDS_FACE. 
2. Attach your favorite manifold to the TopoDS_Face (I’d personally only use 
NormalToMeshProjectionManifold) using cell->set_all_manifold_ids(manifold_id) 
(Notice the use of set_all_manifold_ids, and **not** set_manifold_id: you want 
all children of these objects to belong to the same manifold, in particular you 
want all faces that are internal, i.e., that are shared between two obects with 
the same manifold id, to inherit the same manifold id). 
3. Go one codimension lower: in your case, curves (for 3d meshes, surfaces). 
Follow the same rule as above, setting all_manifold_ids on faces that you know 
should follow a known codimension one geometry (a known TopoDS_EDGE, or 
TopoDS_WIRE). In this case I’d only use ArcLengthProjectionManifold objects, 
attached to the wires that identify your geometry.

Doing things in this order guarantees that internal edges get the correct 
manifold id, which, in your case, is not happening.

Best,
Luca.



> On 1 Apr 2020, at 17:32, Bruno Blais  wrote:
> 
> So my investigation on this issue continues...
> It seems the core of my issues stem from using a non-trivial mesh with more 
> than one cell as a starting point.
> if I start with a 2 cell spacedim=2 dim=3 mesh such as:
> 
> 
> 
> 
> 
> 
> Which is identical to the second step of the adaptation. I get an aberrant 
> result (see attached out_wrong.mp4).
> 
> 
> 
> I also get issues when I try to apply the manifold on a starting msh made 
> with GMSH that has a similar topology such as in the following image:
> 
> 
> 
> 
> 
> In this case, I get an error about a point not being on the manifold (when it 
> clearly is)...
> 
> So I must admit I am at loss here. Right now the code I use to set the 
> manifold is :
> 
> for (const auto  : tria.active_cell_iterators())
>   {
> 
> cell->set_manifold_id(1);
> 
> for (const auto  : cell->face_iterators())
> 
>   {
> 
>   face->set_manifold_id(2);
> 
>   }
> 
>   }
> 
> If anybody has any experience with using the OpenCascade manifold using 
> initial meshes that have more than 1 cell, I would greatly appreciate advices 
> :)!.
> 
> Best
> 
> Bruno
> 
> 
> 
> 
> On Wednesday, April 1, 2020 at 9:38:27 AM UTC-4, Bruno Blais wrote:
> 
> So it seems part of my problem relate to the topology I was trying to adapt 
> OR to the staring complex mesh I was using.
> I have managed to make it work starting with a trivial mesh and using a more 
> simple IGES CAD. I will work on complexifying it from there. I will try to 
> keep the information here as it seems that feature has not been used 
> extensively :)!.
> What I have working now is this animation.
> 
> On Tuesday, March 31, 2020 at 9:19:23 AM UTC-4, Bruno Blais wrote:
> Dear all,
> I hope you are well in these uncertain times.
> 
> I have been trying to use the OpenCASCADE library to set-up my manifolds. The 
> goal of this is to be able to do CFD simulations in complex geometry, but use 
> an initial very coarse mesh made with GMSH and use the dynamics mesh 
> adaptation to adapt the mesh while having the adaptation be 

[deal.II] Difficulty setting manifold with Opencascade using gmsh mesh + salome STEP (or IGES) - Step-54

[Bruno Blais]

Dear all,
I hope you are well in these uncertain times.

I have been trying to use the OpenCASCADE library to set-up my manifolds. 
The goal of this is to be able to do CFD simulations in complex geometry, 
but use an initial very coarse mesh made with GMSH and use the dynamics 
mesh adaptation to adapt the mesh while having the adaptation be CAD-aware. 
The geometry is complex, so fixing the manifold analytically appears to be 
a bad idea.
I started from step-54, which compiles and works well on my machines.
However, I am unable to "make it work" with my own generated Salome STEP 
and gmsh MSH files.
My CAD is in metric, so I adjusted the scaling to use a scaling of 1. My 
mesh also contains numerous wires, so I disabled that aspect and only 
focused on refining with the faces which I label by looping over them.
However, no matter what occurs I get an error thrown:



An error occurred in line <114> of file 
 
in function
dealii::Point 
dealii::OpenCASCADE::NormalProjectionManifold::project_to_manifold(const dealii::ArrayView >&, const dealii::Point&) const [with int 
dim = 2; int spacedim = 3]
The violated condition was: 
closest_point(sh, surrounding_points[i], tolerance) 
.distance(surrounding_points[i]) < std::max(tolerance * 
surrounding_points[i].norm(), tolerance)
Additional information: 
The point [ 0.1 0 0 ] is not on the manifold.

Stacktrace:
---
#0  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::OpenCASCADE::NormalProjectionManifold<2, 
3>::project_to_manifold(dealii::ArrayView const, 
dealii::MemorySpace::Host> const&, dealii::Point<3, double> const&) const
#1  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::FlatManifold<2, 
3>::get_new_point(dealii::ArrayView const, 
dealii::MemorySpace::Host> const&, dealii::ArrayView const&) const
#2  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::Manifold<2, 
3>::get_new_point_on_line(dealii::TriaIterator > const&) const
#3  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
#4  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
#5  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::TriaAccessor<1, 2, 3>::center(bool, bool) const
#6  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::Triangulation<2, 3>::DistortedCellList 
dealii::internal::TriangulationImplementation::Implementation::execute_refinement<3>(dealii::Triangulation<2,
 
3>&, bool)
#7  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::Triangulation<2, 3>::execute_refinement()
#8  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::Triangulation<2, 3>::execute_coarsening_and_refinement()
#9  /home/blab/soft/dealii/dealii/build/lib/libdeal_II.g.so.9.2.0-pre: 
dealii::Triangulation<2, 3>::refine_global(unsigned int)
#10  ./step-54: Step54::TriangulationOnCAD::refine_mesh()
#11  ./step-54: Step54::TriangulationOnCAD::run()
#12  ./step-54: main


I have tried some elements :
- Playing with the tolerance of the CAD
- Writing my CAD from the OpenCASCADE STEP object and then opening it again 
in Salome to see that it is interpreted correctly +

However, I always seem to face an issue where my CAD does not appear to 
generate a valid manifold. I have joined my code as well as the salome HDF 
file and gmsh GEO file used to generate the CAD and meshes respectively.
I must admit I am quite at lost here as to why both don't coincide. Any 
advices on where to proceed from there?
Best :)!
Bruno

-- 
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.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/fdc5899b-53d1-41ca-8efe-72742a52aed6%40googlegroups.com.


topology_adaptation.tar.gz
Description: Binary data