From: Roy Stogner <roystgnr@ic...>  20060626 21:20:49

On Mon, 26 Jun 2006, David Knezevic wrote: > I was wondering what the preferred method is for imposing Dirichlet > BCs on a domain with triangles (TRI6s in my case), My preferred method for setting Dirichlet BCs is by actually integrating the penalized terms over each boundary side. Trying to set BCs node by node on LAGRANGE elements is faster, but trying to do it node by node on any other element is wrong. > since we have boundary nodes that are not on a "boundary side" of an > element. Every boundary node should be on a boundary side of some element  by what other definition is it a boundary node? Even if you're adding a penalty matrix node by node, you don't have to add the penalty on every element containing each boundary node, you just have to add it on some element containing each boundary node. If a boundary node is part of three elements but it's only on a boundary edge of two of them, that's fine  it doesn't matter if the final system contains 2*penalty on the matrix diagonal and rhs or if it contains 3*penalty; as long as the same value is on each, you'll still get the right boundary value to within O(1/penalty) error. > Currently I'm looping over the elements and testing the location of > each node to see if a BC should be imposed using the penalty method, > i.e. I use the following code within the element loop: > > > for(unsigned int n=0; n<elem>n_nodes(); n++) > { > const Point& p = elem>point(n); > if(fabs(p(1)x2_max) < TOL) > { > Kuu(n,n) += penalty; // set u = 0 on this boundary > Kvv(n,n) += penalty; > Fv(n) += penalty; // set v = 1 on this boundary > } > } > > > I think this should work (?), but it's giving me weird behaviour > where the BC value seems to depend on the number of times the > particular node is looped over. Could you be more specific? What BC values are you seeing along y==x2_max? > I was wondering if there's a better way of doing this, and if there's > anything wrong with the approach I described above. The approach above will break on nonLAGRANGE elements and on any elements which use first order bases on a second order geometric element. It assumes that every node at y==x2_max is a boundary node, and it leaves the natural boundary conditions on every other boundary. None of those statements is necessarily a problem, though, so I'm not sure what the bug could be.  Roy 