From: Roy Stogner <roystgnr@ic...>  20110110 19:59:18

On Mon, 20 Dec 2010, robert wrote: > I am trying to build a model which solves the same equation as in ex9. > I would like to specify the boundary conditions in a way that the flux > out of the system on the left side equals the flux coming in on the > right side. > I am trying to use something like: > > > 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(); > > const std::vector<Real>& JxW_face = fe_face>get_JxW(); > > const std::vector<Point >& qface_point = fe_face>get_xyz(); > > fe_face>reinit(elem, side); > > for (unsigned int qp=0; qp<qface.n_points(); qp++) > { > const Real xf = qface_point[qp](0); > const Real yf = qface_point[qp](1); > const Real zf = qface_point[qp](2); > > > > const Real value = neumann_value(xf, yf, zf); > > for (unsigned int i=0; i<phi_face.size(); i++) > Fe(i) += JxW_face[qp]*value*phi_face[i][qp]; > } > > } > but I am not sure how to determine neumann value for the outgoing flux on the other side. Not sure whether this eventually got answered or figured out or not, but basically you want to: Take the xyz value at each quadrature point, Add to it the translation vector between its side and the other side, Call MeshBase::point_locator()() on the result to get the corresponding other_side element. Use FEInterface::inverse_map to get the master element point on the other side, Call fe_other_side>reinit to get flux values there.  Roy 