On Thu, Jul 07, 2005 at 04:30:35PM +0200, Karin Kraft wrote:
Hello!
We have made some progress. We have made a simple implementation
of the
general boundary conditions. As a start we have assumed homogeneous
dirichlet with a constant penalty factor and this works reasonably
well.
Excellent!
Some issues we have discovered:
1. We need to compute the local edge numbers. At the moment we do
this
by comparing to all the edges in the cell like so:
// Iterate over all edges in the boundary
for (EdgeIterator edge(boundary); !edge.end(); ++edge)
{
...
for(int i = 0; i < cell.noEdges(); i++)
{
if(cell.edgeID(i) == edge->id())
{
segment = i;
}
}
However, this is inefficient. Perhaps the mesh should be extended to
store this information?
Matt Knepley is working on a new parallel mesh component for FEniCS
and the plan is that we will wrap this in DOLFIN (same as with PETSc)
once it is ready, so maybe we should not invest too much time in
extending the mesh data structures with new data.
I suggest that we use your temporary implementation (looping over
edges locally to find the local edge number), but create it as a
function in Cell (put it close to the functions for computing edge
and
face alignments). Something like
uint Cell::edgeNumber(const Edge& edge) const;
We will also need
uint Cell::faceNumber(const Face& edge) const;
By putting this as a member function in Cell, you can use the local
data structures directly instead of iterators (the array ce) which
should be a little faster.
2. Our eval looks something like this and we want to avoid the
conditional evaluation. Any ideas?
void eval(real block[], const AffineMap& map, uint boundary) const
{
// Compute geometry tensors
real G0_ = map.scale*1e6;
// Compute element tensor
if( boundary == 0 ) {
block[0] = 0.0;
block[1] = 0.0;
block[2] = 0.0;
block[3] = 0.0;
block[4] = 3.333333333333318e-01*G0_;
block[5] = 3.333333333333318e-01*G0_;
block[6] = 0.0;
block[7] = 3.333333333333318e-01*G0_;
block[8] = 3.333333333333318e-01*G0_;
}
Regards,
Karin and Johan
I managed to avoid conditionals when building the dofmap()
function by
creating arrays. Maybe there is some solution where the segment
(boundary) argument is used to index an array to get the values, but
it doesn't look like that would be easy to do here.
One solution would be to create three different eval functions
(eval_segment_0, eval_segment_1. eval_segment_2) and then have an
array of function pointers to these functions and use the argument
segment (boundary) to index into this array, get the correct function
and then call it. I'm not sure how this compares to having a
switch or
a series of if statements inside eval.
It bugs me that we first have to search for the edge number, and then
when we have the number go through a switch statement. Are we doing
something wrong? Maybe there is a shortcut that avoids both the
searching and the switch?
/Anders
_______________________________________________
DOLFIN-dev mailing list
DOLFIN-dev@xxxxxxxxxx
http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev
