From: Roy Stogner <roystgnr@ic...>  20060426 19:28:06

On Wed, 26 Apr 2006, David Knezevic wrote: > I was just wondering what exactly happens in libMesh when you set a periodic > boundary condition constraint in order to constrain, say, dof0 and dof1 to > be the same, e.g.: > > DofConstraintRow constraint; > constraint[dof0]=1; > dof_map.add_constraint_row(dof1, constraint); > > Does this simply reset the row in the system matrix corresponding to dof1 to > be [0 0 0 ... 0 1 0 ... 0 1 0 ... 0], where the dof0'th entry is 1 and the > dof1'th entry of the row is 1? Is this all one needs to do to impose > periodic boundary conditions in a finite element code? (of course the RHS > here would need to have a 0 in the dof1'th row as well) No  if you just reset a row like this, you would lose information about the equation that was originally assembled on that row. Periodic BCs are just like adaptivity constraints in this sense: when you constrain dof0 to be a linear multiple of dof1, you are effectively adding dof1's test function to dof0. > My understanding of periodic BCs is that, when imposing periodic BCs, we can > think of two corresponding boundary degrees of freedom on opposite sides of > the domain as the same DOF, and then the support of this DOF's basis > function is split across the boundary of the domain. In this case, in order > to take into account the split support of the basis function, I would think > we would need to do the full matrix assembly, then add row dof1 to row dof0, > then finally set the constraint row as in the above code snippet. Is this in > fact what is required when setting periodic BCs? have a hunch this adding of > matrix rows is unnecessary, clarification would be much appreciated. Adding constraints does require this sort of addition conceptually (otherwise you might be setting up weighted residual equations with weird discontinuous halftestfunctions), but in libMesh we do it element by element. Basically, the user assembles the problem K*x = f, but the constraints tell us that x = C*y, so we instead solve the problem C^T*K*C*y = C^T*f (multiplying on both sides to preserve symmetry). And, just as K is the sum of a bunch of small dense element matrices K_e, so is C, so we can do the multiplications elementbyelement. Take a look at src/base/dof_map_constraints.C for the implementation.  Roy 