Menu
▾
▴
Re: [Libmesh-users] LibMesh applications to phase ordering problems
|
From: John P. <jwp...@gm...> - 2008-04-28 17:50:04
|
On Mon, Apr 28, 2008 at 9:02 AM, Nasser Mohieddin Abukhdeir <nas...@mc...> wrote: > Hi: > Thanks, I ended up figuring that out over the weekend, my background > in C is not serving me very well with C++! > So far I have run into two major difficulties which I think require > me to make some code additions. Note, I am posting this to the dlist in > case I am using LibMesh incorrectly, plus, I plan on eventually > submitting this code to the LibMesh developers if all goes well: > > 1) The AdaptiveTimeSolver implementation is not working well for my > system of PDEs. I think there are two reasons behind this: > a) the method involves solving three nonlinear steps per > overall timestep > which is extremely expensive for my system (on the order of > millions of DOFs) b) my system is "very" nonlinear, so the > provided algorithm is too optimistic and I end up having > to set the max_growth variable to limit the predicted > timestep. > The solution, I think, is to implement a more simple one-step algorithm > based on a discrete PI-control method (see the reference) using the LTE > estimate from the non-linear step. Basically I would end up creating a > new AdaptiveTimeSolver class where you can chose the time-adaptive > algorithm and I would add the one I just mentioned plus maybe another > simple one-step algorithm. Sounds great! The current adaptive timesolver implementation is useful in many of the problems we've tackled but we don't expect it to work perfectly for all problems. Hopefully it will be possible to write your new solver(s) in the existing TimeSolver framework. > 2) Using adaptive meshing is a huge advantage, but I have an issue > finding a good initial mesh. I plan on implementing a simple routine > that repeats the initial time step, re-projecting the initial solution > on the refined mesh. This would get around the problem I have now where > my initial coarse mesh has to be not-so-coarse to resolve my initial > solution (which involves a heaviside function). I'm not sure where this > would fit in to the LibMesh class structure. This is certainly an issue we commonly face. One common solution is to start with a coarse grid and adapt to the initial data; it sounds as though this is what you are already doing. The "exact" error indicator in LibMesh may help you here, if you are not already using it. > ref- "Adaptive Multilevel Solutions of Nonlinear Parabolic PDE Systems," > Jens Lang, Springer 2001, pages 32-35 -J > > > Nasser Mohieddin Abukhdeir > Graduate Student (A.D. Rey Research Group) > McGill University - Department of Chemical Engineering > http://webpages.mcgill.ca/students/nabukh/web/ > http://people.mcgill.ca/alejandro.rey/ > > > > > > John Peterson wrote: > > Hi, > > > > Since your initialization function is a member function of a class you > > can't pass a pointer to it for project_solution the same way you would > > for a "global" function. An easy fix is to use a global function > > instead, but I believe there may be another, more correct way to > > achieve what you are trying to do. > > > > -J > > > > On Wed, Apr 23, 2008 at 3:14 PM, Nasser Mohieddin Abukhdeir > > <nas...@mc...> wrote: > > > >> Thanks for the informative answer, I moved over the DiffSystem and I > >> think I worked out some incorrect math, but now I am back to square one > >> with the initial condition of the transient problem. I added this line > >> to the end of my Nem2DSystem:init_data () call (Nem2DSystem is the name > >> of my class that inherits FEMSystem ): > >> > >> this->project_solution(initial_value, > >> NULL,(this->get_equation_systems()).parameters); > >> > >> where the initial value function is defined similar to the > >> non-DiffSystem example with an initial condition: > >> > >> Number Nem2DSystem::initial_value (const Point& p, > >> const Parameters& parameters, > >> const std::string& sysname, > >> const std::string& unkname) > >> > >> and I get a compile error that I cannot figure out: > >> > >> nem2dsystem.C: In member function \u2018virtual void > >> Nem2DSystem::init_data(): > >> nem2dsystem.C:149: error: no matching function for call to > >> \u2018Nem2DSystem::project_solution(<unresolved overloaded function > >> type>, NULL, Parameters&) > >> /home/grad/nasser/nemfem/libmesh-0.6.2/include/solvers/system.h:165: > >> note: candidates are: void System::project_solution(Number (*)(const > >> Point&, const Parameters&, const std::string&, const std::string&), > >> Gradient (*)(const Point&, const Parameters&, const std::string&, const > >> std::string&), Parameters&) const > >> > >> I guess I have to study C++ templates a little more, but any ideas? > >> > >> Nasser Mohieddin Abukhdeir > >> Graduate Student (A.D. Rey Research Group) > >> McGill University - Department of Chemical Engineering > >> http://webpages.mcgill.ca/students/nabukh/web/ > >> http://people.mcgill.ca/alejandro.rey/ > >> > >> > >> > >> > >> > >> Roy Stogner wrote: > >> > > >> > On Wed, 23 Apr 2008, Nasser Mohieddin Abukhdeir wrote: > >> > > >> >> *** Solving time step 0, time = 1e-07 *** > >> >> ERROR: LASPACK Error: > >> >> in ILUFactor for Mat: > >> >> Factorization produces zero pivot elements. > >> >> [0] src/numerics/laspack_linear_solver.C, line 276, compiled Apr 21 2008 > >> >> at 14:01:19 > >> >> Aborted > >> >> > >> >> Basically I added my expressions for the residual and Jacobian and > >> >> removed all of the boundary condition code because I need to impose > >> >> Neuman BCs (which from what I understand is implicitely taken care of by > >> >> the FEM algorithm). > >> > > >> > Homogeneous natural BCs are implicitly taken care of by the FEM > >> > formulation, but just what the "natural boundary condition" is depends > >> > on your precise formulation; if you've got something complicated, take > >> > a look at what boundary integrals pop out of your integration by parts > >> > to double-check that what you're setting to zero is what you want to > >> > be setting to zero. > >> > > >> >> I initialized t=0 with an acceptable solution to the equations, so I > >> >> just expected some smooth sailing... > >> > > >> > You're solving a transient problem, then? In an all-Neumann problem > >> > with no reaction terms you can often get a singular matrix doing a > >> > steady solve (since the problem itself isn't uniquely defined); but > >> > for a transient problem, the mass matrix coming from the time > >> > integration should prevent anything like that from happening. > >> > > >> >> I am also wondering about DiffSystem, should I just scrap my work with > >> >> ex13 and start with DiffSystem? > >> > > >> > There are tradeoffs. DiffSystem makes it easier to add some new > >> > features (e.g. adaptive time stepping, more complex nonlinear > >> > equations) to your code, but it also means you'll have to make > >> > more changes to your code to keep up with the library. > >> > > >> > When you construct a DiffSystem subclass, you read: > >> > "Warning, This code is untested, experimental, or likely to see future > >> > API changes" > >> > > >> > The "untested" and "experimental" adjectives no longer apply, but > >> > "future API changes" certainly does. I'm probably going to alter even > >> > the most basic class methods to make it work better multithreaded. > >> > The changes should be straightforward (methods will take reference > >> > parameter(s) instead of using certain member variables), but it's > >> > never fun when your code suddenly requires repair because a new > >> > library version made backwards incompatible API changes. > >> > --- > >> > Roy > >> > >> ------------------------------------------------------------------------- > >> This SF.net email is sponsored by the 2008 JavaOne(SM) Conference > >> Don't miss this year's exciting event. There's still time to save $100. > >> Use priority code J8TL2D2. > >> http://ad.doubleclick.net/clk;198757673;13503038;p?http://java.sun.com/javaone > >> _______________________________________________ > >> Libmesh-users mailing list > >> Lib...@li... > >> https://lists.sourceforge.net/lists/listinfo/libmesh-users > >> > >> > > ------------------------------------------------------------------------- > This SF.net email is sponsored by the 2008 JavaOne(SM) Conference > Don't miss this year's exciting event. There's still time to save $100. > Use priority code J8TL2D2. > http://ad.doubleclick.net/clk;198757673;13503038;p?http://java.sun.com/javaone > _______________________________________________ > Libmesh-users mailing list > Lib...@li... > https://lists.sourceforge.net/lists/listinfo/libmesh-users > |