Thank you very much for the pointers, Daniel! - Including the templates header did the trick. I no longer need to copy the data over - I was able to implement the penalty method and the DoF constraint method; the subspace projection method seems unambiguously best in terms of accuracy, iterations, and wall clock time. I benchmarked all three approaches against an MS for the pure Neumann problem. - Re: spectrum--no, I only need the max and min eigenvalues to get at the conditioning of the system. However, this isn't crucial as I can just measure how the number of conjugate gradient iterations scale as a proxy.
Cheers, Corbin On Friday, May 19, 2023 at 1:47:42 PM UTC-4 d.arnd...@gmail.com wrote: > Corbin, > > - Using ARPACK might indeed be an option. Do you really need to know the > whole spectrum, though? > - Yes, you would just call AffineConstraints::add_line() to add constrain > a given DoF to zero. Note that it's a no-op if the DoF is already > constrained. This doesn't change the size of the linear system. > - You will need to include deal.II/lac/affine_constraints.templates.h, > deal.II/lac/affine_constraints.h only contains the declarations but not the > definitions. > > Best, > Daniel > > On Fri, May 19, 2023 at 12:04 PM Corbin Foucart <corbin....@gmail.com> > wrote: > >> I'm working on a similar pure Neumann problem with a rank-1 deficiency >> (pressure known only up to a constant). I've implemented the subspace >> projection method Konrad mentioned, with good results. However, I'm >> interested in comparing it to the single DoF constraint method, as well as >> an alternative penalization method alluded to in the paper by Bochev. I had >> a few questions: >> >> 1. For a Sparse or ChunkSparse matrix object, is there a simple way >> of obtaining the eigenvalues? I know I can rebuild deal.ii with ARPACK, >> but >> this seems like overkill >> 2. Exactly how would I use an AffineConstraints object to "pin" the >> degree of freedom? Is it a matter of a single call to >> AffineConstraints::add_line()? It's unclear to me from the documentation >> whether this would add an additional "row" to the linear system, or >> replace >> an existing one (I'd be looking to do the second) >> 3. I've implemented the subspace projection by wrapping the >> ChunkSparse matrix class and overriding the vmult() function; however, >> at >> the moment, I'm forced to copy an original ChunkSparse matrix object into >> the wrapped object because the >> AffineConstraints::distribute_local_to_global() function doesn't >> recognize >> the wrapped matrix --I receive a linker error >> >> libHDG_experimental.so: error: undefined reference to 'void >> dealii::AffineConstraints::distribute_local_to_global<HDG::LinearSystem::ChunkSparseWrapper, >> >> dealii::Vector >(dealii::FullMatrix const&, dealii::Vector const&, >> std::vector<unsigned int, std::allocator > const&, >> HDG::LinearSystem::ChunkSparseWrapper&, dealii::Vector&, bool, >> std::integral_constant<bool, false>) const' >> >> I find that surprising, given that my wrapper class inherits >> publicly from ChunkSparseMatrix. >> >> Any guidance would be appreciated! >> Corbin >> >> On Monday, December 28, 2020 at 1:23:19 PM UTC-5 Wolfgang Bangerth wrote: >> >>> On 12/26/20 3:06 AM, Konrad Simon wrote: >>> > What one can also do is just constrain one DoF to a specific value >>> (this would >>> > also remove rigid motion in elasticity). But think about your solution >>> > variable: If it is in the Sobolev space H^1 then point evaluations may >>> not be >>> > defined for dimension larger than 2. Similarly if, for example, the >>> pressure >>> > in a mechanical or fluid problem is often just in L^2. Point >>> evaluations do >>> > not make sense there at all. >>> >>> Right, this is the correct approach: Constrain a single degree of >>> freedom to >>> zero (or any other value you choose) and solve the problem. Then you can >>> subtract the mean value of the solution *after* solving the linear >>> system. >>> (See VectorTools::subtract_mean_value and >>> VectorTools::compute_mean_value.) >>> >>> If you're uncomfortable with the ill-posedness of taking a point value, >>> you >>> can also take the mean value along a small segment of the boundary >>> (step-1) or >>> a small part of the domain. But in practice, this is not necessary and >>> Konrad's solution is what everyone seems to be doing. >>> >>> Best >>> W. >>> >>> -- >>> ------------------------------------------------------------------------ >>> Wolfgang Bangerth email: bang...@colostate.edu >>> www: http://www.math.colostate.edu/~bangerth/ >>> >>> -- >> The deal.II project is located at http://www.dealii.org/ >> For mailing list/forum options, see >> https://groups.google.com/d/forum/dealii?hl=en >> --- >> You received this message because you are subscribed to the Google Groups >> "deal.II User Group" group. >> To unsubscribe from this group and stop receiving emails from it, send an >> email to dealii+un...@googlegroups.com. >> To view this discussion on the web visit >> https://groups.google.com/d/msgid/dealii/ae39fbcc-c76b-4c03-bfbb-93d19ab2e3d7n%40googlegroups.com >> >> <https://groups.google.com/d/msgid/dealii/ae39fbcc-c76b-4c03-bfbb-93d19ab2e3d7n%40googlegroups.com?utm_medium=email&utm_source=footer> >> . >> > -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/5639d536-91e8-44e1-a3dc-aad23c3129e5n%40googlegroups.com.
