The periodic boundary condition for an RVE analysis

[Farzin Mozafari]

 

Dear all,

 

I am a new user of Deal ii. I am going to solve a multiscale problem in solid mechanics using Deal ii. To this end, I am going to model a simple cube with a linear elastic material under periodic boundary condition as the RVE of the material containing its microstructural features. Then by applying simple shear and uniaxial strain-controlled boundary conditions, I am going to determine the associated components of the averaged stress and elasticity tensor.

 

The way that usually people use is to involve an H matrix (as explained in the attached file). In the literature, this method is called the generalized periodic boundary condition. Indeed, this method makes us able to directly apply shear and uniaxial loading through the H matrix components (i.e. Displacement gradient matrix). I am not sure how can I implement this method in deal II.

 

I would be grateful if anyone gives me a piece of advice to carry it out properly.

 

Regards

Farzin

 

 

[Daniel Arndt]

Farzin,

it appears what you want to do is almost what DoFTools::make_periodicity_constraints(https://www.dealii.org/current/doxygen/deal.II/namespaceDoFTools.html#a929249499b1e5624728d212e90a8e037) gives you.

You "only" have to apply the additional displacement +H(\chi_+-\chi_-). You can check for the support points of the constraints using DoFTools::map_dofs_to_support_points(https://www.dealii.org/current/doxygen/deal.II/namespaceDoFTools.html#a5514e4f59ea659f63953d62ca429eaff) assuming that you are using a nodal finite element. With this information, you should be able to calculate the displacement for each constrained dof pair

and add it using AffineConstraints::set_inhomogeneity (https://www.dealii.org/current/doxygen/deal.II/classAffineConstraints.html#a370b2b485a03fe2bb82c1098a6bdfc4c).

Best,

Daniel

--
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+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/7D60CE96-E78D-45FA-A5C7-C51AAD2E3918%40hxcore.ol.

[Jean-Paul Pelteret]

Hi Farzin,

To complement what Daniel already said, I noted in this previous post an approach whereby you solve for the microscopic fluctuation field directly (i.e. using homogeneous periodic constraints) and then correct for the superimposed linear deformation field as a post-processing step. Although I’d previously approached the problem in the same general manner that you’ve described (specifically, introducing a Lagrangian periodicity frame to compute the inhomogeneity for the affine constraints), a colleague of mine found success using the simpler method of first solving for the fluctuation field. 

I hope that this (well, the information in my prior post) helps a bit. 

Best,

Jean-Paul

To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/CAOYDWbLbQmfrgPJxsoYc8k%2B54vF6rugJ%2B4hHuA61mfDv9HKA3A%40mail.gmail.com.

[Farzin Mozafari]

Dear Daniel,

Thank you for your reply.

I tried to follow your solution, but I am encountering an error. To what follows, I will post my code and error as well. I will be thankful if you give me a piece of advice to get rid of this error. 
As I mentioned, my domain is a linear elastic cube that I am going to impose PBCs on opposite facec.

My code is:

template <int dim>

void CubeVector<dim>::meshing()

{

const Point<dim> LL(0, 0, 0);

const Point<dim> RH(0.1, 0.1, 0.1);

std::vector<unsigned int> repetition(dim); 

repetition[0] = 1;

repetition[1] = 1;

repetition[2] = 1;

GridGenerator::subdivided_hyper_rectangle(triangulation, repetition, LL, RH);

std::string mesh_name = "Rectangle-mesh";

std::ofstream outputmeshfile(mesh_name);

GridOut gridout;

gridout.write_eps(triangulation, outputmeshfile);

std::cout << " Mesh file output to "

<< mesh_name

<< std::endl;

// Boundary indicator // Assign name to every face. This indicators will be used for defining periodic faces.

double Lenght_x = RH[0] - LL[0];

double Lenght_y = RH[1] - LL[1];

double Lenght_z = RH[2] - LL[2];

for (auto& face : triangulation.active_face_iterators()){

if (face->at_boundary())

{

if (face->center()[0] == 0) { face->set_boundary_id(10);}

if (face->center()[0] == Lenght_x) { face->set_boundary_id(11);}

if (face->center()[1] == 0) { face->set_boundary_id(20);}

if (face->center()[1] == Lenght_y) { face->set_boundary_id(21);}

if (face->center()[2] == 0) { face->set_boundary_id(30);}

if (face->center()[2] == Lenght_z) { face->set_boundary_id(31);}

}

        }

std::vector<GridTools::PeriodicFacePair<

GridGenerator::Triangulation<dim>::cell_iterator>>

periodicity_vector;

GridTools::collect_periodic_faces(triangulation,

10,

11,

0,

periodicity_vector);

}

And the error is :

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc: In member function ‘void CubeVector<dim>::meshing()’:

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:179:19: error: ‘Triangulation’ is not a member of ‘dealii::GridGenerator’

  179 |    GridGenerator::Triangulation<dim>::cell_iterator>>

      |                   ^~~~~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:179:19: note: suggested alternatives:

In file included from /usr/local/include/deal.II/grid/tria_description.h:23,

                 from /usr/local/include/deal.II/grid/tria.h:28,

                 from /mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:2:

/usr/local/include/deal.II/grid/cell_id.h:37:7: note:   ‘dealii::Triangulation’

   37 | class Triangulation;

      |       ^~~~~~~~~~~~~

In file included from /usr/local/include/deal.II/distributed/tria.h:27,

                 from /usr/local/include/deal.II/fe/fe_tools.h:30,

                 from /usr/local/include/deal.II/fe/fe_system.h:28,

                 from /mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:11:

/usr/local/include/deal.II/distributed/tria_base.h:302:9: note:   ‘dealii::parallel::Triangulation’

  302 |   using Triangulation DEAL_II_DEPRECATED = TriangulationBase<dim, spacedim>;

      |         ^~~~~~~~~~~~~

In file included from /usr/local/include/deal.II/fe/fe_tools.h:30,

                 from /usr/local/include/deal.II/fe/fe_system.h:28,

                 from /mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:11:

/usr/local/include/deal.II/distributed/tria.h:1403:11: note:   ‘dealii::parallel::distributed::Triangulation’

 1403 |     class Triangulation

      |           ^~~~~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:179:36: error: template argument 1 is invalid

  179 |    GridGenerator::Triangulation<dim>::cell_iterator>>

      |                                    ^

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:179:39: error: template argument 1 is invalid

  179 |    GridGenerator::Triangulation<dim>::cell_iterator>>

      |                                       ^~~~~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:179:39: error: template argument 2 is invalid

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:179:52: error: expected unqualified-id before ‘>’ token

  179 |    GridGenerator::Triangulation<dim>::cell_iterator>>

      |                                                    ^~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:182:14: error: ‘collect_periodic_faces’ is not a member of ‘dealii::GridTools’

  182 |   GridTools::collect_periodic_faces(triangulation,

      |              ^~~~~~~~~~~~~~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:186:4: error: ‘periodicity_vector’ was not declared in this scope

  186 |    periodicity_vector);

      |    ^~~~~~~~~~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc: At global scope:

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:190:70: error: template argument 1 is invalid

  190 |  std::vector<GridTools::PeriodicFacePair<Triangulation::cell_iterator>> periodicity_vector;

      |                                                                      ^~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:190:73: error: template argument 1 is invalid

  190 |  std::vector<GridTools::PeriodicFacePair<Triangulation::cell_iterator>> periodicity_vector;

      |                                                                         ^~~~~~~~~~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:190:73: error: template argument 2 is invalid

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:192:39: error: expected constructor, destructor, or type conversion before ‘(’ token

  192 |      GridTools::collect_periodic_faces(dof_handler,10,11,0,matched_pairs);

      |                                       ^

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc: In instantiation of ‘void CubeVector<dim>::setup_system() [with int dim = 3]’:

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:377:2:   required from ‘void CubeVector<dim>::run() [with int dim = 3]’

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:387:17:   required from here

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:198:21: warning: unused variable ‘total_dofs’ [-Wunused-variable]

  198 |  const unsigned int total_dofs = dof_handler.n_dofs();

      |                     ^~~~~~~~~~

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc: In instantiation of ‘void CubeVector<dim>::assemble() [with int dim = 3]’:

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:378:2:   required from ‘void CubeVector<dim>::run() [with int dim = 3]’

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:387:17:   required from here

/mnt/d/research/deal_ii_projects/3d_elasticity/step-1.cc:289:9: warning: unused variable ‘x’ [-Wunused-variable]

  289 |  double x = fe_face_values.quadrature_point(face_q)[0];

      |         ^

make[3]: *** [CMakeFiles/step-1.dir/build.make:63: CMakeFiles/step-1.dir/step-1.cc.o] Error 1

make[2]: *** [CMakeFiles/Makefile2:272: CMakeFiles/step-1.dir/all] Error 2

make[1]: *** [CMakeFiles/Makefile2:252: CMakeFiles/run.dir/rule] Error 2

make: *** [Makefile:196: run] Error 2

Thank you so much in advance!

[Daniel Arndt]

Farzin,

std::vector<GridTools::PeriodicFacePair<

GridGenerator::Triangulation<dim>::cell_iterator>>

periodicity_vector;

should be

std::vector<GridTools::PeriodicFacePair<

Triangulation<dim>::cell_iterator>>

periodicity_vector;

In general, I would recommend looking into step-45 (https://www.dealii.org/current/doxygen/deal.II/step_45.html) for an example using periodic boundary conditions.

Best,

Daniel

To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/d2248f66-6c74-495c-8a79-0b69120212e7n%40googlegroups.com.

[Farzin Mozafari]

Dear Daniel,

Thank you for your reply.

As per your recommendation, I followed step-45 for applying PBC. As far as I realized, since I am not using parallel meshing, the second and third steps are not necessary for my case. Thus, I just used steps 1 and 4 as follows:

        DoFTools::make_periodicity_constraints(dof_handler,10,11,0,ComponentMask(),constraints,1);

But, I am getting the following error:

error: no matching function for call to ‘make_periodicity_constraints(dealii::DoFHandler<3, 3>&, int, int, int, dealii::ComponentMask, dealii::AffineConstraints<double>&, int)’

  188 |         DoFTools::make_periodicity_constraints(dof_handler,10,11,0,ComponentMask(),constraints,1);

      |         ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~^~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

Are there any further steps that need to be taken?

Regards,

Farzin 

[Daniel Arndt]

Farzin,

https://www.dealii.org/current/doxygen/deal.II/namespaceDoFTools.html#a03324e6e060c6a55191b2c60ad871eab shows that the signature of the function you want to call is

void DoFTools::make_periodicity_constraints(

  const DoFHandlerType &   dof_handler,

  const types::boundary_id   b_id1,

  const types::boundary_id   b_id2,

  const unsigned int   direction,

  AffineConstraints< number > &   constraints,

  const ComponentMask &   component_mask = ComponentMask(),

  const number   periodicity_factor = 1.)

In particular, the AffineConstraints object comes before the ComponentMask object (which can actually be omitted in your case).

Best,

Daniel

To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/54bee2b2-f613-400d-a520-cb9432c14462n%40googlegroups.com.