On Thu, Dec 5, 2013 at 3:22 PM, Roy Stogner <[email protected]>wrote:
> > > > >> My experience with HIERARCHIC has been unfavorable, with >> convergence issues of iterative solvers. I have not seen any such >> problems with SZABAB, and BERNSTEIN, so am favoring those. >> > > Strange. From a brief glance at SZABAB it looked like they were > practically just HIERARCHIC with different scaling factors, the kind > of thing even Jacobi preconditioning would fix. > > I remember triying all sorts of preconditioning, and nothing helped. My understanding is that the HIERARCHIC polynomials get highly ill-conditioned as the p-order increases. Infact Karniadakis, in his book on spectral hp-elements, has a plot of condition number for HIERARCHIC, LAGRANGE, and LEGENDRE polynomials for increasing p, and shows that HIERARCHIC are the first to go through the roof, followed by LAGRANGE. I am interested in implementing a C0 continuous version fo LEGENDRE, which modifies the first two basis from a constant and linear to a combination of two linears, and maintains the bubble functions. The orthogonality of these polynomials can be quite useful, especially in getting a diagonal mass matrix. > > I managed to get the higher order SZABAB working with three >> modifications: >> >> -- modified FEInterface::extra_hanging_dofs() to be true for SZABAB. >> >> -- FEInterface::compute_constraints() does not do anything for SZABAB. >> Added a case SZABAB: for that. >> >> -- fe_szabab.C needs FE<2, SZABAB>::compute_constraints. >> >> With this, h-refinement of higher order SZABAB works without problems. >> > > Great! Send us a PR? > > > I am not sure if this is the issue plaguing hp-refinement, but I >> will get to it soon. >> > > No, sadly you ran into another bug on the way to the hp bug. > --- > Roy aahhh... ok... I will continue pushing forward... -Manav ------------------------------------------------------------------------------ Sponsored by Intel(R) XDK Develop, test and display web and hybrid apps with a single code base. Download it for free now! http://pubads.g.doubleclick.net/gampad/clk?id=111408631&iu=/4140/ostg.clktrk _______________________________________________ Libmesh-users mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/libmesh-users
