Hello,
I have been having some issues when trying to run a "toy" code (a simple
Steady Navier-Stokes solver) on a case I used to run.
The boundary group Id 0, which is manually set in my GMSH script
(vonKarman.geo), is detected perfectly as both a BC and a manifold when I
export to GMSH 2 and I get the expected results.
However, if I export the same GMSH file to V4 the same code does not run
anymore. I guess there is a change in behavior I did not notice in how
dealii interprets the new format? What am I doing wrong?
The error I get is (when running in debug mode):
--------------------------------------------------------
An error occurred in line <10387> of file </home/blaisb/codes/dealii/dealii/
source/grid/tria.cc> in function
void dealii::Triangulation<dim, spacedim>::
set_all_manifold_ids_on_boundary(dealii::types::boundary_id, dealii::types::
manifold_id) [with int dim = 2; int spacedim = 2; dealii::types::boundary_id
= unsigned int; dealii::types::manifold_id = unsigned int]
The violated condition was:
boundary_found
Additional information:
The given boundary_id 0 is not defined in this Triangulation!
Stacktrace:
-----------
#0 /home/blaisb/codes/dealii/build/lib/libdeal_II.g.so.9.1.0-pre:
dealii::Triangulation<2, 2>::set_all_manifold_ids_on_boundary(unsigned int,
unsigned int)
#1 ./schursteadynavierstokes: SchurNavierStokesSolver<2>::runVonKarman()
#2 ./schursteadynavierstokes: main
--------------------------------------------------------
However, BC ID 0 as a physical group is defined explicitely in the script.
I have joined the .h and .cc of my code ( I am sorry it is ugly as hell, it
was just made for a quick test).
The GMSH file is as follow (also in the tar archive):
y0=0;
y1=1;
x0 =0.;
x1 =10;
xc=1;
yc=0.5;
r=0.025;
lo = 2.0e-1;
lc = 1.0e-1;
Point(0) = {x0, y0, 0, lo};
Point(1) = {x0, y1, 0, lo};
Point(2) = {x1, y0, 0, lo};
Point(3) = {x1, y1, 0, lo};
Point(4) = {xc, yc, 0, lc};
Point(5) = {xc-r, yc, 0, lc};
Point(6) = {xc, yc+r, 0, lc};
Point(7) = {xc+r, yc, 0, lc};
Point(8) = {xc, yc-r, 0, lc};
// Contour Box
Line(1) = {0,1};
Line(2) = {1,3};
Line(3) = {3,2};
Line(4) = {2,0};
// Inner Circle
Circle(5)={5,4,6};
Circle(6)={6,4,7};
Circle(7)={7,4,8};
Circle(8)={8,4,5};
Line Loop(1) = {1,2,3,4};
Line Loop(2) = {5,6,7,8};
// Creates the physical entities
Physical Line(0)={5,6,7,8};
Physical Line(1)={1};
Physical Line(2)={2,4};
Physical Line(3)={3};
Plane Surface(1) = {1,2};
Physical Surface(2) = {1};
Recombine Surface{1};
The DEALII lines regarding the manifold are :
grid_in.read_msh(input_file);
Point<dim,double> circleCenter(1.0,0.5);
static const SphericalManifold<dim> boundary(circleCenter);
triangulation.set_all_manifold_ids_on_boundary(0,0);
triangulation.set_manifold (0, boundary);
And the lines linked to the BC are :
emplate <int dim>
void SchurNavierStokesSolver<dim>::setup_dofs ()
{
system_matrix.clear();
pressure_mass_matrix.clear();
dof_handler.distribute_dofs(fe);
std::vector<unsigned int> block_component(dim+1, 0);
block_component[dim] = 1;
DoFRenumbering::component_wise (dof_handler, block_component);
dofs_per_block.resize (2);
DoFTools::count_dofs_per_block (dof_handler, dofs_per_block,
block_component);
unsigned int dof_u = dofs_per_block[0];
unsigned int dof_p = dofs_per_block[1];
FEValuesExtractors::Vector velocities(0);
{
nonzero_constraints.clear();
DoFTools::make_hanging_node_constraints(dof_handler,
nonzero_constraints);
VectorTools::interpolate_boundary_values(dof_handler,
0,
ZeroFunction<dim>(dim+1),
nonzero_constraints,
fe.component_mask(velocities
));
if (simulationCase_==TaylorCouette)
{
VectorTools::interpolate_boundary_values(dof_handler,
1,
RotatingWall<dim>(),
nonzero_constraints,
fe.component_mask(
velocities));
}
if (simulationCase_==BackwardStep)
{
VectorTools::interpolate_boundary_values(dof_handler,
1,
PoiseuilleInlet<dim>(),
nonzero_constraints,
fe.component_mask(
velocities));
}
if (simulationCase_==vonKarman)
{
std::set<types::boundary_id> no_normal_flux_boundaries;
no_normal_flux_boundaries.insert (2);
VectorTools::compute_no_normal_flux_constraints (dof_handler, 0,
no_normal_flux_boundaries,
nonzero_constraints
);
}
if (simulationCase_==vonKarman)
{
VectorTools::interpolate_boundary_values(dof_handler,
1,
ConstantXInlet<dim>(),
nonzero_constraints,
fe.component_mask(
velocities));
}
}
nonzero_constraints.close();
{
zero_constraints.clear();
DoFTools::make_hanging_node_constraints(dof_handler, zero_constraints
);
VectorTools::interpolate_boundary_values(dof_handler,
0,
ZeroFunction<dim>(dim+1),
zero_constraints,
fe.component_mask(velocities
));
if (simulationCase_==TaylorCouette || simulationCase_==BackwardStep ||
simulationCase_==vonKarman)
{
VectorTools::interpolate_boundary_values(dof_handler,
1,
ZeroFunction<dim>(dim+1),
zero_constraints,
fe.component_mask(velocities
));
}
if (simulationCase_==vonKarman)
{
std::set<types::boundary_id> no_normal_flux_boundaries;
no_normal_flux_boundaries.insert (2);
VectorTools::compute_no_normal_flux_constraints (dof_handler, 0,
no_normal_flux_boundaries,
zero_constraints
);
}
}
zero_constraints.close();
std::cout << " Number of active cells: "
<< triangulation.n_active_cells()
<< std::endl
<< " Number of degrees of freedom: "
<< dof_handler.n_dofs()
<< " (" << dof_u << '+' << dof_p << ')'
<< std::endl;
}
--
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.
src.tar.gz
Description: Binary data
