data out shows segmentation fault for calculated stress when using more than 8 ranks
59 views
Skip to first unread message
Omkar Kolekar
May 17, 2025, 11:12:53 PM5/17/25
to deal.II User Group
Respected Sir (/Ma'am),
I thank you Prof. W. Bangerth and team for developing Dealii library.
I am trying to solve a two material homogenization problem in linear elastic regime using Dealii. I have a RVE with a material1 and spherical inclusion of material2. Using dealii I was able to calculate the solution, strain and stresses. I am running the code in parallel on 8 ranks.
The problem I encounter is that the code fails when I assign the values from stress to a dealii::vector for visualization. When I run the code on 4 ranks it works fine but for 8 ranks it shows segmentation faults. I used the debugger to check what is causing the problem and I have the following output...
Checking size of s_11: 207425 I am rank 6
Checking local_index: 0 I am rank 6
Checking if there is a problem with the Stress vector:
0.00873843 0.000175943 -6.5817e-05
0.000175943 -7.20853e-05 0.000405372
-6.5817e-05 0.000405372 0.000508171
Thread 1 "simulate_physic" received signal SIGSEGV, Segmentation fault.
Thread 1 "simulate_physic" received signal SIGSEGV, Segmentation fault.
dealii::ElasticProblem<3>::output_results (this=0x7fffffff99b0, seed=4, sample_id=3, target_path=0x7fffffff9800, direction=0, shear_direction=35) at /home/ok28rume/Thesis/thesis/code/modules/linear_elasticity/src/write_results.cpp:346
346 s_12[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][1];
Catchpoint 1 (throw)
dealii::ElasticProblem<3>::output_results (this=0x7fffffff99b0, seed=4, sample_id=3, target_path=0x7fffffff9800, direction=0, shear_direction=35) at /home/ok28rume/Thesis/thesis/code/modules/linear_elasticity/src/write_results.cpp:350
350 s_22[cell->active_cell_index()] = Stress[cell->active_cell_index()][1][1];
#0 dealii::ElasticProblem<3>::output_results (this=0x7fffffff99b0, seed=4, sample_id=3, target_path=0x7fffffff9800, direction=0, shear_direction=35) at /home/ok28rume/Thesis/thesis/code/modules/linear_elasticity/src/write_results.cpp:346
cell = @0x7ffffffed5d0: {<dealii::TriaIterator<dealii::CellAccessor<3, 3> >> = {<dealii::TriaRawIterator<dealii::CellAccessor<3, 3> >> = {accessor = {<dealii::TriaAccessor<3, 3, 3>> = {<dealii::TriaAccessorBase<3, 3, 3>> = {static structure_dimension = 3, present_level = 8, present_index = 65927Catchpoint 1 (throw)
Clearly the Stress (which is a vector of tensors) is a problem here. However the same code runs for 4 ranks but gives errors for 8 ranks.
Checking size of s_11: 207425 I am rank 6
Checking local_index: 0 I am rank 6
Checking if there is a problem with the Stress vector:
0.00873843 0.000175943 -6.5817e-05
0.000175943 -7.20853e-05 0.000405372
-6.5817e-05 0.000405372 0.000508171
Thread 1 "simulate_physic" received signal SIGSEGV, Segmentation fault.
Thread 1 "simulate_physic" received signal SIGSEGV, Segmentation fault.
dealii::ElasticProblem<3>::output_results (this=0x7fffffff99b0, seed=4, sample_id=3, target_path=0x7fffffff9800, direction=0, shear_direction=35) at /home/ok28rume/Thesis/thesis/code/modules/linear_elasticity/src/write_results.cpp:346
346 s_12[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][1];
Catchpoint 1 (throw)
dealii::ElasticProblem<3>::output_results (this=0x7fffffff99b0, seed=4, sample_id=3, target_path=0x7fffffff9800, direction=0, shear_direction=35) at /home/ok28rume/Thesis/thesis/code/modules/linear_elasticity/src/write_results.cpp:350
350 s_22[cell->active_cell_index()] = Stress[cell->active_cell_index()][1][1];
#0 dealii::ElasticProblem<3>::output_results (this=0x7fffffff99b0, seed=4, sample_id=3, target_path=0x7fffffff9800, direction=0, shear_direction=35) at /home/ok28rume/Thesis/thesis/code/modules/linear_elasticity/src/write_results.cpp:346
cell = @0x7ffffffed5d0: {<dealii::TriaIterator<dealii::CellAccessor<3, 3> >> = {<dealii::TriaRawIterator<dealii::CellAccessor<3, 3> >> = {accessor = {<dealii::TriaAccessor<3, 3, 3>> = {<dealii::TriaAccessorBase<3, 3, 3>> = {static structure_dimension = 3, present_level = 8, present_index = 65927Catchpoint 1 (throw)
Clearly the Stress (which is a vector of tensors) is a problem here. However the same code runs for 4 ranks but gives errors for 8 ranks.
I tried to search on the internet but I did not find anything related to this, hence writing here on the forum.
Below is the code as well.
I thank you all in advance for your kind help, suggestions, time and support.
Warm Regards,
Omkar
Post processing code
#include "simulate_physics.h"
template class dealii::ElasticProblem<3>;
namespace dealii {
template <int dim>
void ElasticProblem<dim>::postprocessing() {
pcout<<"***************************************************** LINEAR ELASTICITY MODULE :: PostProcessing Started *********************************************"<<"\n"<< std::endl;
SymmetricTensor<4, dim> local_stiffness_tensor;
Tensor<2, dim> temp_Epsilon; // Temporary strain tensor
Tensor<2, dim> temp_stress; // Temporary average stress tensor
// set the component extractor for the vector values problems
FEValuesExtractors::Vector u_fe;
static const unsigned int first_u_component = 0;
u_fe = first_u_component; //selection starts from the first component.
// Get the number of locally owned cells
const unsigned int n_locally_owned_cells = triangulation.n_locally_owned_active_cells();
// Resizing the vectors of stress, strain, and the displacement gradient.
//TODO: This size may be too much and cause a problem may be??!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
// Resizing the vectors of stress, strain, and the displacement gradient.
F.resize(n_locally_owned_cells);
Epsilon.resize(n_locally_owned_cells);
Stress.resize(n_locally_owned_cells);
//std::cout<<"PostProcessing :: The size of F is " << F.size() << " the size of Epsilon is " << Epsilon.size() << " and the size of Stess is " << Stress.size() << " I am rank " << get_rank() << std::endl;
// Defining the quadrature formula
const QGauss<dim> quadrature_formula(fe.degree + 1);
const std::vector<double> &q_weights = quadrature_formula.get_weights();
// FEValues setup
FEValues<dim> fe_values(fe,
quadrature_formula,
update_values | update_gradients |
update_quadrature_points | update_JxW_values);
// [gradient of all displacements] for each quadrature point in the cell
const unsigned int n_q_points = quadrature_formula.size();
// Create storage for gradients at quadrature points
std::vector<Tensor<2, dim>> grad_u_q(n_q_points);
// Storage for deformation gradient tensor
Tensor<2, dim> grad_u;
barF = 0.0;
barEpsilon = 0.0;
barStress = 0.0;
// Ensure all ranks have completed the above operations
//MPI_Barrier(mpi_communicator);
unsigned int local_cell_index = 0;
for (const auto &cell : dof_handler.active_cell_iterators()){
//std::cout<<"Postprocessing :: I entered the for loop to iterated over active cells its the first for loop. I am rank " << Utilities::MPI::this_mpi_process(mpi_communicator) <<" \n";
if(cell->is_locally_owned()){
fe_values.reinit(cell);
grad_u = 0.0;
//std::cout<<"Postprocessing :: This cell is locally owed by rank " << Utilities::MPI::this_mpi_process(mpi_communicator) <<" \n";
// Check material ID for the current cell
unsigned int material_id = cell->material_id();
// Formulate the stiffness tensor.
if (material_id == 2 ) {
local_stiffness_tensor = get_stress_strain_tensor(lambda2, mu2);
} else {
local_stiffness_tensor = get_stress_strain_tensor(lambda1, mu1);
}
// Get gradients * solution on each quadrature point
fe_values[u_fe].get_function_gradients(solution, grad_u_q);
// Sum contributions of quadrature points
double sum_JxW = 0.0; // The Jacobian is not used anywhere in the code further, hence commented.
for (unsigned int q=0; q < n_q_points; q++){
grad_u += grad_u_q[q] * q_weights[q]; // This is not used anywhere in the code further.
sum_JxW += 1 * fe_values.JxW(q); // The Jacobian is not used anywhere in the code further, hence commented.
}
// Compute Deformation gradient
F[local_cell_index] = Physics::Elasticity::Kinematics::F(grad_u);
barF += cell->measure() * F[local_cell_index];
// Compute linear strain tensor
Epsilon[local_cell_index] = Physics::Elasticity::Kinematics::epsilon(grad_u);
barEpsilon += cell->measure() * Epsilon[local_cell_index];
temp_Epsilon = Epsilon[local_cell_index];
// Compute Stresses.
TensorAccessors::contract<2,4,2,dim>(temp_stress, local_stiffness_tensor, temp_Epsilon);
Stress[local_cell_index] = temp_stress;
int indexi_counter = 0;
int indexj_counter;
//###################################### debugger #########################################################//
for (unsigned int i = 0; i < dim; ++i){
indexj_counter = 0;
for (unsigned int j = 0; j < dim; ++j){
//std::cout << "Postprocessing :: Stress tensor check for local index: " << local_cell_index << " and i = "<< i << " and j = "<< j <<" is " << temp_stress[i][j] << " I am rank " << get_rank() << std::endl;
double value = Stress[local_cell_index][i][j];
++indexj_counter;
}
if (indexj_counter != dim){
std::cout<< "Postprocessing :: WARNING! :: the index did not run for dim times, i is " << indexi_counter << " and the j index is " << indexj_counter << " and the rank is " << get_rank() << std::endl;
}
++indexi_counter;
}
if (indexi_counter != dim){
std::cout<< "Postprocessing :: WARNING! :: the index did not run for dim times, i is " << indexi_counter << " and the j index is " << indexj_counter << " and the rank is " << get_rank() << std::endl;
}
//#########################################################################################################//
barStress += cell->measure() * temp_stress;
local_cell_index++;
}
}
//std::cout<<"Postprocessing :: I exited the for loop. I am rank " << Utilities::MPI::this_mpi_process(mpi_communicator) <<" \n";
// Ensure all ranks have completed the above operations
//MPI_Barrier(mpi_communicator);
//std::cout<<"Postprocessing :: Setting up the global vectors. I am rank " << Utilities::MPI::this_mpi_process(mpi_communicator) <<" \n";
// Compute global averages using MPI_Reduce
Tensor<2, dim> global_barF;
Tensor<2, dim> global_barEpsilon;
Tensor<2, dim> global_barStress;
// Ensure all ranks have completed the above operations
//MPI_Barrier(mpi_communicator);
//std::cout<<"Postprocessing :: Setting the global vectors done. " << Utilities::MPI::this_mpi_process(mpi_communicator) <<" \n";
// Reduce the local bar variables to the global bar variables
//TODO: This is not working as expected. The global variables are not being updated.
/*MPI_Reduce(&barF, &global_barF, sizeof(Tensor<2, dim>), MPI_BYTE, MPI_SUM, 0, mpi_communicator);
MPI_Reduce(&barEpsilon, &global_barEpsilon, sizeof(Tensor<2, dim>), MPI_BYTE, MPI_SUM, 0, mpi_communicator);
MPI_Reduce(&barStress, &global_barStress, sizeof(Tensor<2, dim>), MPI_BYTE, MPI_SUM, 0, mpi_communicator);*/
/*global_barF /= Utilities::MPI::n_mpi_processes(mpi_communicator);
global_barEpsilon /= Utilities::MPI::n_mpi_processes(mpi_communicator);
global_barStress /= Utilities::MPI::n_mpi_processes(mpi_communicator);*/
global_barF = Utilities::MPI::sum(barF, mpi_communicator);
global_barEpsilon = Utilities::MPI::sum(barEpsilon, mpi_communicator);
global_barStress = Utilities::MPI::sum(barStress, mpi_communicator);
//MPI_Barrier(mpi_communicator);
// Copy the global_bar variables back to the local bar variables on the main rank
if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0) {
barF = global_barF;
barEpsilon = global_barEpsilon;
barStress = global_barStress;
}
//MPI_Barrier(mpi_communicator);
//std::cout<<"Postprocessing :: Ended postprocessing by all ranks cause this is after a barrier. I am rank " << Utilities::MPI::this_mpi_process(mpi_communicator) <<" \n";
pcout<<"**************************************************** LINEAR ELASTICITY MODULE :: PostProcessing Completed ********************************************"<<"\n"<< std::endl<< std::endl;
}
template <int dim>
SymmetricTensor<4, dim> ElasticProblem<dim>::get_stress_strain_tensor(const double lambda,
const double mu){
SymmetricTensor<4, dim> tmp;
for (unsigned int i = 0; i < dim; ++i){
for (unsigned int j = 0; j < dim; ++j){
for (unsigned int k = 0; k < dim; ++k){
for (unsigned int l = 0; l < dim; ++l){
tmp[i][j][k][l] = (((i == k) && (j == l) ? mu : 0.0) +
((i == l) && (j == k) ? mu : 0.0) +
((i == j) && (k == l) ? lambda : 0.0));
}
}
}
}
return tmp;
}
}
Write results code (Segmentation fault at s_12[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][1];)
void ElasticProblem<dim>::output_results(int seed, int sample_id, const std::string* target_path, int direction, int shear_direction) {
const unsigned int n_locally_owned_cells = triangulation.n_locally_owned_active_cells();
if(F.size() == 0){
std::cerr << "LINEAR_ELASTICITY :: No results to output. F is empty on rank " << get_rank() << ". " << std::endl;
return;
}
if(Epsilon.size() == 0){
std::cerr << "LINEAR_ELASTICITY :: No results to output. Epsilon is empty on rank " << get_rank() << ". " << std::endl;
return;
}
if(Stress.size() == 0){
std::cerr << "LINEAR_ELASTICITY :: No results to output. Stress is empty on rank " << get_rank() << ". " << std::endl;
return;
}
//#################################################################### Assertions ####################################################################//
AssertThrow(Stress.size() == Epsilon.size(), dealii::ExcMessage("Stress and Epsilon must be the same size!"));
AssertThrow(Stress.size() == F.size(), dealii::ExcMessage("Stress and F must be the same size!"));
AssertThrow(F.size() == Epsilon.size(), dealii::ExcMessage("F and Epsilon must be the same size!"));
std::cout << "Stress size: " << Stress.size() << " F size: " << F.size() << " Epsilon size: "<< Epsilon.size()<< " rank: " << get_rank() << std::endl;
//####################################################################################################################################################//
std::string direction_str;
if(shear_direction == 35){
direction_str = "uniaxial_";
switch (direction) {
case 0: direction_str += "xx_";
break;
case 1: direction_str += "yy_";
break;
case 2: direction_str += "zz_";
break;
default: std::cerr << "Invalid direction for uniaxial loading." << std::endl;
return;
}
} else{
direction_str = "shear_";
switch (direction) {
case 0: direction_str += "x";
break;
case 1: direction_str += "y";
break;
case 2: direction_str += "z";
break;
default: std::cerr << "Invalid direction for shear loading." << std::endl;
return;
}
switch (shear_direction) {
case 0: direction_str += "x_";
break;
case 1: direction_str += "y_";
break;
case 2: direction_str += "z_";
break;
default: std::cerr << "Invalid shear direction." << std::endl;
return;
}
}
//MPI_Barrier(mpi_communicator);
move_mesh();
//MPI_Barrier(mpi_communicator);
DataOut<dim> data_out;
data_out.attach_dof_handler(dof_handler);
// Step 1: Create a vector to hold the material IDs
dealii::Vector<float> material_ids(triangulation.n_active_cells());
//solution.update_ghost_values(); //Updating ghost values as it may cause a problem.
//dealii::Vector<float> material_ids(n_locally_owned_cells);
//############################## DEBUGGER OUTPUTS ##############################
// Step 2: Loop over all cells and store their material IDs
unsigned int cell_index = 0;
for (const auto &cell : triangulation.active_cell_iterators()) {
if (cell->is_locally_owned()){
material_ids[cell_index] = cell->material_id();
++cell_index;
}
}
// Step 3: Attach the material_id as cell data
data_out.add_data_vector(material_ids, "material_id");
std::vector<std::string> solution_names;
switch (dim){
case 1:
solution_names.emplace_back("displacement");
break;
case 2:
solution_names.emplace_back("x_displacement");
solution_names.emplace_back("y_displacement");
break;
case 3:
solution_names.emplace_back("x_displacement");
solution_names.emplace_back("y_displacement");
solution_names.emplace_back("z_displacement");
break;
default:
DEAL_II_NOT_IMPLEMENTED();
}
data_out.add_data_vector(solution, solution_names);
//##################### Adding subdomains ########################//
Vector<float> subdomain(triangulation.n_active_cells());
for (unsigned int i = 0; i < subdomain.size(); ++i)
subdomain(i) = triangulation.locally_owned_subdomain();
data_out.add_data_vector(subdomain, "subdomain");
//######################## Kinematics ##########################//
//***************** Gradient of Displacement *******************//
dealii::Vector<double> F_11(triangulation.n_active_cells());
dealii::Vector<double> F_12(triangulation.n_active_cells());
dealii::Vector<double> F_13(triangulation.n_active_cells());
dealii::Vector<double> F_21(triangulation.n_active_cells());
dealii::Vector<double> F_22(triangulation.n_active_cells());
dealii::Vector<double> F_23(triangulation.n_active_cells());
dealii::Vector<double> F_31(triangulation.n_active_cells());
dealii::Vector<double> F_32(triangulation.n_active_cells());
dealii::Vector<double> F_33(triangulation.n_active_cells());
//************************** Strain ****************************//
dealii::Vector<double> e_11(triangulation.n_active_cells());
dealii::Vector<double> e_12(triangulation.n_active_cells());
dealii::Vector<double> e_13(triangulation.n_active_cells());
dealii::Vector<double> e_21(triangulation.n_active_cells());
dealii::Vector<double> e_22(triangulation.n_active_cells());
dealii::Vector<double> e_23(triangulation.n_active_cells());
dealii::Vector<double> e_31(triangulation.n_active_cells());
dealii::Vector<double> e_32(triangulation.n_active_cells());
dealii::Vector<double> e_33(triangulation.n_active_cells());
//************************** Stress ****************************//
dealii::Vector<double> s_11(triangulation.n_active_cells());
dealii::Vector<double> s_12(triangulation.n_active_cells());
dealii::Vector<double> s_13(triangulation.n_active_cells());
dealii::Vector<double> s_21(triangulation.n_active_cells());
dealii::Vector<double> s_22(triangulation.n_active_cells());
dealii::Vector<double> s_23(triangulation.n_active_cells());
dealii::Vector<double> s_31(triangulation.n_active_cells());
dealii::Vector<double> s_32(triangulation.n_active_cells());
dealii::Vector<double> s_33(triangulation.n_active_cells());
unsigned int local_index = 0;
for (const auto &cell : triangulation.active_cell_iterators()){
if (cell->is_locally_owned() || cell->is_ghost()) {
/*############################################## DEBUGGING ##########################################################*/
std::cout << "Checking size of s_11: " << s_11.size() << " I am rank "<< get_rank() << std::endl;
std::cout << "Checking local_index: " << local_index << " I am rank "<< get_rank() << std::endl;
std::cout << "Checking if there is a problem with the Stress vector: \n";
/*for (unsigned int i = 0; i < 3; ++i){
for (unsigned int j = 0; j < 3; ++j){
std::cout << Stress[local_index][i][j] << " ";
}
std::cout << std::endl;
}*/
/*#################################################################################################################*/
//******************** Gradient of Displacement (F) **********************//
F_11[cell->active_cell_index()] = F[cell->active_cell_index()][0][0];
F_12[cell->active_cell_index()] = F[cell->active_cell_index()][0][1];
F_13[cell->active_cell_index()] = F[cell->active_cell_index()][0][2];
F_21[cell->active_cell_index()] = F[cell->active_cell_index()][1][0];
F_22[cell->active_cell_index()] = F[cell->active_cell_index()][1][1];
F_23[cell->active_cell_index()] = F[cell->active_cell_index()][1][2];
F_31[cell->active_cell_index()] = F[cell->active_cell_index()][2][0];
F_32[cell->active_cell_index()] = F[cell->active_cell_index()][2][1];
F_33[cell->active_cell_index()] = F[cell->active_cell_index()][2][2];
//************************** Strain (epsilon) ****************************//
e_11[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][0][0];
e_12[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][0][1];
e_13[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][0][2];
e_21[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][1][0];
e_22[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][1][1];
e_23[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][1][2];
e_31[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][2][0];
e_32[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][2][1];
e_33[cell->active_cell_index()] = Epsilon[cell->active_cell_index()][2][2];
//*************************** Stress (sigma) *****************************//
s_11[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][0];
s_12[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][1];
s_13[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][2];
s_21[cell->active_cell_index()] = Stress[cell->active_cell_index()][1][0];
s_22[cell->active_cell_index()] = Stress[cell->active_cell_index()][1][1];
s_23[cell->active_cell_index()] = Stress[cell->active_cell_index()][1][2];
s_31[cell->active_cell_index()] = Stress[cell->active_cell_index()][2][0];
s_32[cell->active_cell_index()] = Stress[cell->active_cell_index()][2][1];
s_33[cell->active_cell_index()] = Stress[cell->active_cell_index()][2][2];
}
}
//std::cout << " LINEAR ELASTICITY MODULE :: Completed for and if loop at line 257. I am rank "<< get_rank() << std::endl; //DEBUGGING
//***** Gradient of Displacement (F) *****//
data_out.add_data_vector(F_11, "F11", dealii::DataOut<dim>::type_cell_data);
//std::cout << " LINEAR ELASTICITY MODULE :: Adding vector F11 to data out on line 261. I am rank "<< get_rank() << std::endl; //DEBUGGING
data_out.add_data_vector(F_12, "F12", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_13, "F13", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_21, "F21", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_22, "F22", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_23, "F23", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_31, "F31", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_32, "F32", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(F_33, "F33", dealii::DataOut<dim>::type_cell_data);
//***** Strain (epsilon) *****//
data_out.add_data_vector(e_11, "e11", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_12, "e12", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_13, "e13", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_21, "e21", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_22, "e22", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_23, "e23", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_31, "e31", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_32, "e32", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(e_33, "e33", dealii::DataOut<dim>::type_cell_data);
//***** Stress (sigma) *****//
data_out.add_data_vector(s_11, "stress_11", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_12, "stress_12", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_13, "stress_13", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_21, "stress_21", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_22, "stress_22", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_23, "stress_23", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_31, "stress_31", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_32, "stress_32", dealii::DataOut<dim>::type_cell_data);
data_out.add_data_vector(s_33, "stress_33", dealii::DataOut<dim>::type_cell_data);
//std::cout << " LINEAR ELASTICITY MODULE :: Building patches at line 286. I am rank "<< get_rank() << std::endl; //DEBUGGING
solution.update_ghost_values(); //Updating ghost values as it may cause a problem.
data_out.build_patches();
//std::cout << " LINEAR ELASTICITY MODULE :: Completed build patches at line 289. I am rank "<< get_rank() << std::endl; //DEBUGGING
//std::string filename = *target_path + "/data/ensemble_"+std::to_string(seed)+"/linear_elasticity/solution_"+std::to_string(sample_id)+".vtk";
//ABOVE WAS USED FOR A NON PARALLEL CASE REFER BELOW FOR PARALLEL CASE.
//std::string filename = *target_path + "/data/ensemble_"+std::to_string(seed)+"/linear_elasticity/Sample_"+std::to_string(sample_id)+direction_str;
std::string filename = *target_path + "/data/ensemble_"+std::to_string(seed)+"/linear_elasticity/"; //THIS IS FOR DEBUGGING
std::string output_name = "Sample_"+std::to_string(sample_id)+"_solution_"+direction_str;
data_out.write_vtu_with_pvtu_record(filename, output_name, 0, mpi_communicator, 2 , 0); //FOR PARALLEL CASE.
//std::ofstream output(filename);
//data_out.write_vtk(output);
//MPI_Barrier(mpi_communicator);
//pcout << " LINEAR ELASTICITY MODULE :: write results line 304. I am rank "<< get_rank() << std::endl;
unmove_mesh();
//MPI_Barrier(mpi_communicator);
if (Utilities::MPI::this_mpi_process(mpi_communicator) == 0) {
save_stress_strain(seed, sample_id, target_path);
}
}
Wolfgang Bangerth
May 19, 2025, 1:02:40 PM5/19/25
to dea...@googlegroups.com
Omkar:
since you are already in the debugger, and know that the problem is in a
statement of the form
s_12[cell->active_cell_index()] = Stress[cell->active_cell_index()][0][1];
the next step would be to investigate whether accessing individual elements of
the arrays in question here is valid. A reasonable suspicion would be that,
for example, s_12 has size 200 and cell->active_cell_index() returns 210, to
give just one example.
I have not tried to look into this too carefully, but see from your code that
you do
Stress.resize(n_locally_owned_cells);
But the number of locally owned cells can be less than
cell->active_cell_index(), and so I am not surprised that you get a
segmentation fault. The statement needs to read
Stress.resize(n_active_cells);
Best
Wolfgang
On 5/16/25 17:39, Omkar Kolekar wrote:
> *** Caution: EXTERNAL Sender ***
> The deal.II project is located at http://www.dealii.org/ <https://
> nam10.safelinks.protection.outlook.com/?
> url=http%3A%2F%2Fwww.dealii.org%2F&data=05%7C02%7CWolfgang.Bangerth%40colostate.edu%7C97f4df0270554d5a38a408dd95b9e1c7%7Cafb58802ff7a4bb1ab21367ff2ecfc8b%7C0%7C0%7C638831347808044058%7CUnknown%7CTWFpbGZsb3d8eyJFbXB0eU1hcGkiOnRydWUsIlYiOiIwLjAuMDAwMCIsIlAiOiJXaW4zMiIsIkFOIjoiTWFpbCIsIldUIjoyfQ%3D%3D%7C0%7C%7C%7C&sdata=L9A%2Bfzn6OJ0qGIpiYGX0sr2QGzKlUQGAeYtVWF%2Bi3mg%3D&reserved=0>
> For mailing list/forum options, see https://groups.google.com/d/forum/dealii?
> hl=en <https://nam10.safelinks.protection.outlook.com/?
> url=https%3A%2F%2Fgroups.google.com%2Fd%2Fforum%2Fdealii%3Fhl%3Den&data=05%7C02%7CWolfgang.Bangerth%40colostate.edu%7C97f4df0270554d5a38a408dd95b9e1c7%7Cafb58802ff7a4bb1ab21367ff2ecfc8b%7C0%7C0%7C638831347808065920%7CUnknown%7CTWFpbGZsb3d8eyJFbXB0eU1hcGkiOnRydWUsIlYiOiIwLjAuMDAwMCIsIlAiOiJXaW4zMiIsIkFOIjoiTWFpbCIsIldUIjoyfQ%3D%3D%7C0%7C%7C%7C&sdata=VYZbiQowwgK%2B8xKtxsFwVSDmsVjRI1XWT3wR5g9dG0U%3D&reserved=0>
> ---
> 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+un...@googlegroups.com
> <mailto:dealii+un...@googlegroups.com>.
> To view this discussion visit https://groups.google.com/d/msgid/dealii/
> ef334937-cf1c-4ae8-93af-4d5752c53800n%40googlegroups.com <https://
> nam10.safelinks.protection.outlook.com/?
> url=https%3A%2F%2Fgroups.google.com%2Fd%2Fmsgid%2Fdealii%2Fef334937-
> cf1c-4ae8-93af-4d5752c53800n%2540googlegroups.com%3Futm_medium%3Demail%26utm_source%3Dfooter&data=05%7C02%7CWolfgang.Bangerth%40colostate.edu%7C97f4df0270554d5a38a408dd95b9e1c7%7Cafb58802ff7a4bb1ab21367ff2ecfc8b%7C0%7C0%7C638831347808082188%7CUnknown%7CTWFpbGZsb3d8eyJFbXB0eU1hcGkiOnRydWUsIlYiOiIwLjAuMDAwMCIsIlAiOiJXaW4zMiIsIkFOIjoiTWFpbCIsIldUIjoyfQ%3D%3D%7C0%7C%7C%7C&sdata=c2dK0UoOSO8D8ETSPIvx7WvT3bgtX9uwkzV1ezH4JfQ%3D&reserved=0>.
Omkar Kolekar
May 20, 2025, 4:18:38 PM5/20/25
to deal.II User Group
Respected Professor Bangerth,
Thank you very much for you help and support.
Stress.resize(n_active_cells);
The above line worked for me.
The above line worked for me.
Warm Regards,
Omkar Kolekar
Reply all
Reply to author
Forward
0 new messages