## libmesh-users

 [Libmesh-users] Is this correct for applying periodic boundary condition? From: Shengli Xu - 2006-02-13 11:19:38 Attachments: Message as HTML ```Dear Buddies, The following code consulted Michael and Ben's article. Uniform meshes, apply periodic boundary on x=3D0 and x=3D1.0 in x direction= . 1, question about Michael's code. ... if (elem->contains_point(point2)) =09=09=09{ =09=09=09 =09=09=09 //elem contains the opposite point =09=09=09 DofMap::DofConstraintRow constraint; ``````````````````````??? Is it correct? =09=09 =09=09=09 dof_map.dof_indices (elem, dof_indices); =09=09=09 =09=09=09 std::vector point2_vec(1); =09=09=09 point2_vec[0] =3D point2; ... for (int i1 =3D 0; i1 < phi.size(); i1++) =09=09=09 { =09=09=09 const std::vector& qface_point =3D fe->get_xyz(); =20 ```````````??? qface_point is not use, what is it for? =09=09=09 if ( std::abs(phi[i1][0]) > func_tol )=20 constraint[dof_indices[i1]] =3D phi[i1][0]; =09=09=09 } 2, The following is my code. But I can not get the correct periodic result on u in the direction x. void apply_periodic_bc(EquationSystems &es, const std::string &system_name)= { assert(system_name=3D=3D"stokes_flow"); PerfLog perf_log("Periodic bc. Assembly",false); const Mesh &mesh=3Des.get_mesh(); LinearImplicitSystem &system=3Des.get_system ("stokes_flow"); DofMap &dof_map=3Dsystem.get_dof_map(); const double min_coord[2]=3D{0.,0.}; const double max_coord[2]=3D{1.,1.}; const double pos_tol=3D1e-10; // x direction periodic for u for(unsigned int n=3D0;npos_tol)){ // get node on the left boundary for(unsigned int n1=3D0;n1
 [Libmesh-users] Re: Is this correct for applying periodic boundary condition? From: Wout Ruijter - 2006-02-27 11:13:41 Attachments: elastic_2D.C ```ShengliXu First of all, you're trying to solve something that doesn't make sense (to me at least), you try to compress a square in x-direction and at the same time demand the left and right bound to have the same x-displacement. You can't get a result that shows the effects of your periodic code because the penalty is 1e10 and therefore the DofConstraints (of order 1) are overshadowed. Always print out the convergence history of your solver, if you have trouble converging usually something like that is happening. If you want the bounds to have the same shape but allow them to move closer together as a whole, you have to add a "ghost" DOF to the system which represents this movement, and add it to all of the dofconstraintrows that you are adding for a certain dimension/dof combination. I needed this before and implemented it by having a loose element added to the system, which will give you a bunch of dofs to use for such purposes (multiple dimension periodicity, for example). Make sure you don't assemble that element though, just set the diagonal terms to nonzero for the dofs you don't use. Secondly, there must be some bug in your code as it only finds one set of opposite dofs, I wrote a periodic BC code in the file that I find works. It might be a matter of taste, but I think the preferred way is still to use element based looping and dof querying as opposed to node based (the smart people on the list might correct me on that ;). Cheers Wout ```
 [Libmesh-users] Re: Is this correct for applying periodic boundary condition? From: Wout Ruijter - 2006-02-28 16:22:38 Attachments: Makefile elastic_2D.C ```ShengliXu // you forgot the all important line: dof_map.constrain_element_matrix_and_vector(Ke, Fe, dof_indices); before system.matrix->add_matrix (Ke, dof_indices); system.rhs->add_vector (Fe, dof_indices); So all the other explanations I came up with for the program not working are far-fetched or nonsense. You can choose the solver on the command line, for example: time /usr/bin/mpirun -np 1 \$(target) -ksp_type gmres -pc_type jacobi >> log= .txt Now, one more time: If you apply a strain in x-direction the bounds are going to move apart, this CONFLICTS with the boundary conditions you just put up (as you're telling them they have the same displacement). As you will see, if you couple all dofs on all opposing sides as you propose you will only get some noise as an answer regardless of the loading. Read my earlier post on how to deal with that (by adding dofs to the system in a dummy element that represent a displacement of a whole side, keeping the sides in the same shape). Secondly, the code you put in to deal with corners does not do anything, as I said, the corners are in the sides and are dealt with automatically (because "couples" is a std::map and keys are unique in those). Cheers W ```