Hello ,
I would like to use finite difference method to calculate the gradient of density in cells.
But I got issue when the program was processed at the part of compute_gradient();
Here is the code below about density and gradient.
//>>>>>>>>>>>>>>>Density>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
template <int dim,typename NumberType>
void Solid<dim,NumberType>::compute_density()
{
typename DoFHandler<dim>::active_cell_iterator cell = dof_handler_ref.begin_active(),
endc = dof_handler_ref.end();
for (; cell != endc; ++cell)
{
// Compute the grid density for each cell,
density[cell->active_cell_index()] = 1.0 / cell->measure();
}
}
//>>>>>>>>>>>>>>>>>>>Gradient>>>>>>>>>>>>>>>>>>>>>>>>
template <int dim,typename NumberType>
void Solid<dim,NumberType>::compute_gradient()
{
typename DoFHandler<dim>::active_cell_iterator cell = dof_handler_ref.begin_active(),
endc = dof_handler_ref.end();
for (; cell != endc; ++cell)
{
for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell; ++face)
{
if (cell->face(face)->at_boundary())
{
// Compute the gradient at each boundary face by differencing the neighboring cell densities
const typename DoFHandler<dim>::cell_iterator neighbor = cell->neighbor(face);
gradient[cell->face(face)->index()] = density[neighbor->active_cell_index()] - density[cell->active_cell_index()];
}
}
}
}
//<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
The error comes out at the part of compute_gradient().
----------------------------------------------------------------------------------------------------------------
Error information:
Thread 1 "cook_membrane" received signal SIGSEGV, Segmentation fault.
0x00007ffff4bfcdaa in dealii::CellAccessor<3, 3>::active_cell_index() const () from /lib/x86_64-linux-gnu/libdeal.ii.so.9.3.2
----------------------------------------------------------------------------------------------------------------
I tried to reduce the numbers of grids and dofs ,but I sill got the same issue.
in the part of void Solid<dim,NumberType>::run(),
I add the code below
template <int dim,typename NumberType>
void Solid<dim,NumberType>::run(){
.//other codes
.//other codes
.//other codes
density.reinit(dof_handler_ref.n_dofs());
gradient.reinit(dof_handler_ref.n_dofs());
//<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
print_vertical_tip_displacement();
std::cout<<"Density claculation"<<std::endl;
compute_density();//sucessfully
std::cout<<"Gradient calculation"<<std::endl;
compute_gradient();//error comes out
}
Could anyone provide any hint or suggestions?
Thanks in advance!
Best regards
Lance