Yes, well this is going to give strange results because you're setting the
boundary ID of internal faces (x=0.5 is the centreline of your geometry).
Before you set any boundary ID's, you should always first check to see that
you're actually on a boundary. The first few tutorials do go through this,
but I've provided you with a worked example to demonstrate the best
practise.
Regards,
J-P
On Sunday, November 13, 2016 at 8:53:01 PM UTC+1, Jaekwang Kim wrote:
>
> Thank you for quick reply!!!
>
> I tried your suggestions but I couldn't resolve the problem.....
>
> The more weird thing is if I commands as..
>
>
> GridGenerator::hyper_cube (triangulation, 0, 1);
>
> triangulation.refine_global ();
>
>
> for (typename Triangulation<dim>::active_cell_iterator
>
> cell=triangulation.begin_active();
>
> cell!=triangulation.end(); ++cell)
>
> {
>
> for (unsigned int f=0;
> f<GeometryInfo<dim>::faces_per_cell; ++f)
>
> {
>
>
>
>
>
> if ( std::abs(cell->face(f)->center()[0]-0.5)< 1e-12 )
>
>
> {
>
>
> cell->face(f)->set_all_boundary_ids(1);
>
>
>
> }
>
>
>
>
>
> }
>
> }
>
>
>
> <https://lh3.googleusercontent.com/-vxXjt_6106w/WCjEWZXkl9I/AAAAAAAAA9M/HtP_S_FZR3wqpaBSHpe7SjgNG9tCnYX5gCLcB/s1600/Screen%2BShot%2B2016-11-13%2Bat%2B1.50.56%2BPM.png>
>
>
>
> 2016년 11월 13일 일요일 오후 1시 35분 3초 UTC-6, Jean-Paul Pelteret 님의 말:
>>
>> Hi Jaekwang,
>>
>> With what you've shown I'm going to guess that your problem might simply
>> be related to numerical / round-off error associated with floating point
>> calculations.
>>
>> if ( cell->face(f)->center()[0]==1 ) // The LHS value is a computed
>>> double, and the RHS value is an integer which will be cast to a double for
>>> comparison
>>
>>
>> Take a look at the mesh-generation part in the introduction to step-2
>> <https://www.dealii.org/developer/doxygen/deal.II/step_2.html#Meshgenerationhttps://www.dealii.org/developer/doxygen/deal.II/step_2.html%23Meshgeneration>
>> for
>> how one best implements these types of comparisons. They're particularly
>> relevant when dealing with the geometry, especially in more complex
>> scenarios. Admittedly the step-2 example is for a circular geometry but the
>> same concept applies:
>>
>> if ( std::abs(cell->face(f)->center()[0]) < 1e-10 ) ...
>>
>> if ( std::abs(cell->face(f)->center()[0] - 1.0) < 1e-10 ) ...
>>
>> You can almost never be certain that any calculated floating point value
>> will be an exact match to an implicitly cast integer value right down to
>> machine precision, which is what will be checked with the equality operator.
>>>
>>>
>> I hope that this fixes your issue.
>>
>> Regards,
>> J-P
>>
>
--
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.
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/grid/tria_iterator.h>
#include <deal.II/grid/grid_generator.h>
#include <iostream>
using namespace dealii;
int main()
{
const int dim = 2;
Triangulation<dim> triangulation;
GridGenerator::hyper_cube (triangulation, 0.0, 1.0);
triangulation.refine_global (1);
std::cout << "Number of active cells: "<< triangulation.n_active_cells()<< std::endl;
for (typename Triangulation<dim>::active_cell_iterator
cell=triangulation.begin_active();
cell!=triangulation.end(); ++cell)
{
for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
{
// You need this line to ensure that this doesn't affect anything other than bondary faces!
if (cell->face(f)->at_boundary() == true)
{
// This is an internal faces! You never want to set this to anything other than its default value.
// The about check should make sure that we never could execute this part of the code.
if ( std::abs(cell->face(f)->center()[0]-0.5)< 1e-12 )
{
cell->face(f)->set_boundary_id(10);
}
// Boundary faces
static const double tol = 1e-12;
if ( std::abs(cell->face(f)->center()[0]-0.0)< tol )
{
cell->face(f)->set_boundary_id(1);
}
if ( std::abs(cell->face(f)->center()[0]-1.0)< tol )
{
cell->face(f)->set_boundary_id(2);
}
if ( std::abs(cell->face(f)->center()[1]-0.0)< tol )
{
cell->face(f)->set_boundary_id(3);
}
if ( std::abs(cell->face(f)->center()[1]-1.0)< tol )
{
cell->face(f)->set_boundary_id(4);
}
}
}
}
int cell_index = 0;
for (typename Triangulation<dim>::active_cell_iterator
cell=triangulation.begin_active();
cell!=triangulation.end(); ++cell, ++cell_index)
{
std::cout << "Cell index: " << cell_index << std::endl;
for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
{
if (cell->face(f)->at_boundary() == true)
std::cout
<< "cell->face(f)->center(): " << cell->face(f)->center()
<< " cell->face(f)->boundary_id(): " << static_cast<int>(cell->face(f)->boundary_id())
<< std::endl;
}
}
return 0;
}
