Brilliant, that would do the job in a very easy way. Ok last question that I would ask related to this topic :-)
I not sure how to find the nodes (DOFs) which live on the boundary (outer surface of my 3D domain). I appreciate if someone can finalize this issue for me. Kindest regards, Omar On Wed, May 11, 2011 at 10:53 AM, Derek Gaston <[email protected]> wrote: > The easiest way to do what you want is by using matrix.zero_rows(). The > way it works is that you fill up a vector of dof indices that correspond to > rows that you want to zero and put a 1 on the diagonal. Then you just call > matrix.zero_rows(dofs_to_zero, 1). The second argument to the function is > the number you want to place on the diagonal (zero by default). > > So.... fill up your matrix like normal, keeping track of rows you want to > wipe out. Close the matrix. Call zero_rows(). Close the matrix again. > > A quick note on this... if you are using Petsc this will modify your > nonzero pattern (ie remove all off-diagonal nonzero entries for these rows). > That's not normally what you want because the next time you come around to > fill the matrix and you do some insertions on that row it will be incredibly > slow (because Petsc will have to reallocate those entries). You can get > beyond this by using something like: > > MatSetOption(static_cast<PetscMatrix<Number> &>(my_matrix).mat(), > MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE); > > That is the correct line for Petsc 3.1. If you are using a different > version there will most likely be a slightly different way to tell it to > keep the nonzero pattern. I recommend setting this each time before you > start filling the matrix up. Petsc seems to like to forget about this > option..... > > Derek > > On May 10, 2011, at 5:10 PM, Omar Al-Abbasi wrote: > > > Thank you all for all the suggestions you made. > > In fact I tried doing two seperate sweeps as Boyce suggested, but it > seems > > that I am doing a mistake which leads to no change in mass matrix. > > > > This how I am doing that after the first sweep for assembling the > matrices: > > . > > . > > . > > for ( ; el != end_el; ++el) > > { > > const Elem* elem = *el; > > dof_map.dof_indices (elem, dof_indices); > > fe->reinit (elem); > > Ke.resize (dof_indices.size(), dof_indices.size()); > > Me.resize (dof_indices.size(), dof_indices.size()); > > > > for (unsigned int side=0; side<elem->n_sides(); side++){ > > if (elem->neighbor(side) == NULL) > > { > > const std::vector<std::vector<Real> >& phi_face = > > fe_face->get_phi(); > > fe_face->reinit(elem, side); > > for (unsigned int i=0; i<phi_face.size(); i++){ > > matrix_B.set(dof_indices[i],dof_indices[i],1.); > > } // end if (elem->neighbor(side) == NULL) > > } // end boundary condition section > > } > > } > > > > Also, what is the best way to loop over the nodes which live on the outer > > surface of the domain? > > Is it just the way I am doing it for (unsigned int side=0; > > side<elem->n_sides(); side++){ > > if (elem->neighbor(side) == NULL) > > and then for (unsigned int i=0; i<phi_face.size(); i++){ .... > > > > Thanks again, > > Omar > > > ------------------------------------------------------------------------------ > > Achieve unprecedented app performance and reliability > > What every C/C++ and Fortran developer should know. > > Learn how Intel has extended the reach of its next-generation tools > > to help boost performance applications - inlcuding clusters. > > http://p.sf.net/sfu/intel-dev2devmay > > _______________________________________________ > > Libmesh-users mailing list > > [email protected] > > https://lists.sourceforge.net/lists/listinfo/libmesh-users > > ------------------------------------------------------------------------------ Achieve unprecedented app performance and reliability What every C/C++ and Fortran developer should know. Learn how Intel has extended the reach of its next-generation tools to help boost performance applications - inlcuding clusters. http://p.sf.net/sfu/intel-dev2devmay _______________________________________________ Libmesh-users mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/libmesh-users
