From: Roy S. <roy...@ic...> - 2004-11-23 00:55:09
|
In the system_projection.C file there's a half-finished implementation of project_vector() which sets up a global L2 projection problem. This is entirely commented out in favor of the local interpolation code above it... but the local interpolation only works for Lagrange elements. Would it be reasonable to split the difference, and make the non-Lagrange projection solve local L2 projection problems on each element? This wouldn't be as efficient as element-specific interpolation code, but it should give the correct results on elements for which perfect interpolation is possible, and should give adequate results for any element. I'm testing this idea myself (it appears to work, and my ex14 works with hierarchic elements now), but I've run into a couple small problems: the NumericVector class has the function: add_vector (const DenseVector< T > &V, const std::vector< unsigned int > &dof_indices) but I couldn't find any equivalent set_vector function. Setting vector elements one by one can't be very efficient if it's a distributed vector. I seem to be having trouble with the new iterators. The following code works: MeshBase::element_iterator elem_it = _mesh.active_local_elements_begin(); const MeshBase::element_iterator elem_end = _mesh.active_local_elements_end(); But if I try to declare these as MeshBase::const_element_iterator type instead (as is done in src/solvers/exact_solution.C), I get the error upon construction: error: conversion from `variant_filter_iterator<Elem*, Predicates::multi_predicate>' to non-scalar type `variant_filter_iterator<Elem* const, Predicates::multi_predicate>' requested If the quadrature rule in project_vector is too low, the element mass matrix may not be SPD and a cholesky solve may fail. At the moment I'm using get_order on the finite element object to figure out what degree of gaussian quadrature is necessary, but this brings up a question: is get_order used to return the highest degree polynomial a finite element uses (important for quadrature) or the highest degree polynomial an element can exactly represent (important for expected convergence rates)? Usually the two are the same, but that may become a problem in the future. --- Roy Stogner |
From: John P. <pet...@cf...> - 2004-11-23 01:54:30
|
Roy Stogner writes: > > > I seem to be having trouble with the new iterators. The following > code works: > > MeshBase::element_iterator elem_it = > _mesh.active_local_elements_begin(); > const MeshBase::element_iterator elem_end = > _mesh.active_local_elements_end(); > > > But if I try to declare these as MeshBase::const_element_iterator type > instead (as is done in src/solvers/exact_solution.C), I get the error > upon construction: > > error: conversion from `variant_filter_iterator<Elem*, > Predicates::multi_predicate>' to non-scalar type > `variant_filter_iterator<Elem* const, Predicates::multi_predicate>' > requested Hi, Only a const Mesh object can construct const_element_iterators, since this requires the calling of a const member function... this is similar to the behavior of STL objects in that foo.begin() returns a const_iterator if foo is a const object. The code in exact_solution.C uses the const_iterator since it has a const reference to the Mesh as private data. Let me know if that clears anything up. -John |
From: Roy S. <roy...@ic...> - 2004-11-23 02:29:45
|
On Mon, 22 Nov 2004, John Peterson wrote: > Only a const Mesh object can construct const_element_iterators, > since this requires the calling of a const member function... You can call const member functions on non-const objects; you just can't call non-const member functions on const objects. Things may be different with functions that are overloaded based on constness, though; I'm not sure how to call the const version of a function which also has a non-const counterpart with the same name and same input arguments. > this is similar to the behavior of STL objects in that foo.begin() > returns a const_iterator if foo is a const object. The code in > exact_solution.C uses the const_iterator since it has a const > reference to the Mesh as private data. Let me know if > that clears anything up. That explains my error; I think it's something that ought to be fixed eventually, though. You can initialize a const pointer to const from a pointer to non-const, and you can initialize an STL const_iterator from a non-const iterator; you ought to be able to do the same with libMesh iterators. --- Roy |
From: John P. <pet...@cf...> - 2004-11-23 02:47:51
|
Roy Stogner writes: > On Mon, 22 Nov 2004, John Peterson wrote: > > > Only a const Mesh object can construct const_element_iterators, > > since this requires the calling of a const member function... > > You can call const member functions on non-const objects; you just > can't call non-const member functions on const objects. Things may be > different with functions that are overloaded based on constness, > though; I'm not sure how to call the const version of a function which > also has a non-const counterpart with the same name and same input > arguments. > > > this is similar to the behavior of STL objects in that foo.begin() > > returns a const_iterator if foo is a const object. The code in > > exact_solution.C uses the const_iterator since it has a const > > reference to the Mesh as private data. Let me know if > > that clears anything up. > > That explains my error; I think it's something that ought to be > fixed eventually, though. You can initialize a const pointer to const > from a pointer to non-const, and you can initialize an STL > const_iterator from a non-const iterator; you ought to be able to do > the same with libMesh iterators. Apparently, not portably and efficiently. Here's a quote from clc++m (search for: construct iterator const_iterator) > Item 26 and Item 27 from Effective STL by Scott Meyers goes into detail > on this problem. There is no portable way of using any casts to convert > from a const_iterator to an iterator. > > Your only option is to create a non const iterator, and advance its > position until it points to the same element as the const_iterator. > > In a vector this will be relatively efficient. In a list, this is an > O(n) operation. Of course, we need to go the other way round, but I imagine you face the same problem. -John |
From: Roy S. <roy...@ic...> - 2004-11-23 03:31:35
|
On Mon, 22 Nov 2004, John Peterson wrote: > Roy Stogner writes: > > Apparently, not portably and efficiently. Here's a quote from clc++m > (search for: construct iterator const_iterator) > > > Item 26 and Item 27 from Effective STL by Scott Meyers goes into detail > > on this problem. There is no portable way of using any casts to convert > > from a const_iterator to an iterator. > Of course, we need to go the other way round, but I imagine you face > the same problem. No, in theory there shouldn't be any problem at all. Converting from a const_iterator to an iterator is a hack that generally means something is wrong with your code; and based on discussions I've read on compiler optimization, the fact that such conversions are even theoretically possible means something is wrong with the C++ standard. Converting from an iterator to a const_iterator, on the other hand, doesn't break any constness assumptions and from a user's perspective should just require you to write "i = j". From a library author's perspective, though, C++ doesn't know that a const_iterator is really an iterator-to-const until you explicitly write conversion operators and such to get it to behave like one. The following code all compiles: vector<int>::iterator i; vector<int>::const_iterator j = i; int const *i; int *j = i; The following does not: MeshBase::element_iterator i; MeshBase::const_element_iterator j = i; --- Roy |