So if I loop over all non-artificial cells instead of locally owned cells to build the dirichlet map, then it seems to work correctly.
Best praveen On Mon, Sep 5, 2016 at 7:37 PM, Praveen C <[email protected]> wrote: > Dear all > I have to create a constraint matrix to apply dirichlet bc by taking the > values from an existing solution (in vector “x" below) on the same mesh. > > for(typename DoFHandler<dim>::active_cell_iterator > cell = dof_handler.begin_active(), > endc = dof_handler.end(); > cell!=endc; ++cell) > if(cell->is_locally_owned()) > { > for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; > ++f) > if (cell->face(f)->at_boundary()) > { > cell->face(f)->get_dof_indices(dof_indices); > for(unsigned int i=0; i<dofs_per_face; ++i) > { > // search if index exists in boundary map > const unsigned int global_i = dof_indices[i]; > std::map<types::global_dof_index,double>::iterator > it = boundary_values_x.find(global_i); > if(it == boundary_values_x.end()) // Did not find > index, so add it > { > boundary_values_x.insert(std: > :pair<types::global_dof_index,double>(global_i,x(global_i))); > } > } > } > } > > add_dirichlet_constraints (boundary_values_x, constraints_x); > constraints_x.close(); > > The function add_dirichlet_constraints is same as here > > https://groups.google.com/forum/#!searchin/dealii/add_ > dirichlet_constraints%7Csort:relevance/dealii/r9C8n1aU3zc/DmtfwWkDAgAJ > > This is a parallel code using Trilinos MPI vectors. I have doubts about > the correctness of this procedure. > > In particular, should I include ghost dofs while creating the boundary > condition map. I am looping over locally owned cells and collecting > whatevar boundary dofs I find into the map. > > Best > praveen > -- 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 [email protected]. For more options, visit https://groups.google.com/d/optout.
