Thank you Derek! I actually looked into MOOSE before my own implementation but then I found this:
https://groups.google.com/forum/#!searchin/moose-users/central$20difference%7Csort:date/moose-users/FJnjdc7u6Xc/juFZpdXMBwAJ Let me follow up the thread on MOOSE and continue there :) Best, Shawn On Tue, May 21, 2019 at 5:37 AM Derek Gaston <fried...@gmail.com> wrote: > Yep - vector_mult, reciprocal, and pointwise_mult are how we do explicit > updating in MOOSE: > > > https://github.com/idaholab/moose/blob/next/framework/src/timeintegrators/ActuallyExplicitEuler.C#L182 > > Derek > > > On Mon, May 20, 2019 at 2:38 PM Yuxiang Wang <yw...@virginia.edu> wrote: > >> Thank you Jed & John! That's extremely helpful. >> >> John - thank you so much for the insight and pointer to reciprocal() and >> pointwise_mult()! As for adding additional NumericVector to the system - >> do >> you happen to know any example or code snippet that did this? What I am >> doing now is creating a new ExplicitSystem and using its solution vector >> as >> the storage, which handles the initialization quite nicely but maybe is an >> overkill. >> >> Best, >> Shawn >> >> On Mon, May 20, 2019 at 7:11 AM Jed Brown <j...@jedbrown.org> wrote: >> >> > John Peterson <jwpeter...@gmail.com> writes: >> > >> > > On Sun, May 19, 2019 at 9:15 PM Yuxiang Wang <yw...@virginia.edu> >> wrote: >> > > >> > >> Dear all, >> > >> >> > >> Sorry for the spam and being new to the field. >> > >> >> > >> I am currently trying to implement an elastodynamics problem with >> > explicit >> > >> method (central difference method, to be specific). I am planning to >> use >> > >> lumped mass matrix (and Rayleigh damping when needed), so the system >> > matrix >> > >> will be simply a diagonal matrix. Solving is therefore trivial - I >> just >> > >> need to do per-element division of the rhs by the diagonal entries to >> > get >> > >> the solution vector. >> > >> >> > >> For this problem, I have the option of just treating it as an >> implicit >> > >> system - fill the system.matrix with only diagonal components and >> call >> > the >> > >> PETSc LU solver. This is very easy thanks to a lot of the libmesh >> > >> infrastructure available. If I do so, should I be concerned about a >> > >> significant slowdown? Or would PETSc be smart enough to realize that >> > this >> > >> is already diagonal matrix and be efficient in solving it? >> > >> >> > > >> > > Yes, I'm pretty sure this will be significantly slower than doing a >> > > "manual" inversion by taking the reciprocal of the matrix diagonal. I >> > don't >> > > know of anything in PETSc's LU solver that will detect this particular >> > > special case. For an explicit code where you want every timestep to >> be as >> > > fast as possible, it will likely be prohibitively slow. >> > >> > If you have a sparse matrix with only the diagonal nonzero, then LU will >> > create no fill and it'll actually be pretty fast (likely hard to measure >> > compared to the cost of evaluating the explicit part), but -pc_type >> > jacobi would also be an exact solver and is more precisely what you >> > need. Note that the default bjacobi/ilu is also an exact solver in this >> > circumstance, as are many other preconditioners. >> > >> > > My other choice would be to create a NumericVector myself to store the >> > >> diagonal system matrix entries, and perform the per-element division. >> > This >> > >> would take more work and will not be using the already well-tested >> > libmesh >> > >> infrastructure, so I am trying to see whether I can get away with >> doing >> > >> this without compromising on the performance. >> > > >> > > >> > > This would probably work best. You can add additional vectors to >> Systems >> > > and assemble into them similarly to the way you assemble the >> right-hand >> > > side vector. Also note that NumericVectors have the reciprocal() API, >> > which >> > > will allow you to quickly compute the inverse, as well as >> > pointwise_mult() >> > > which should allow you to quickly apply it. >> > > >> > > -- >> > > John >> > > >> > > _______________________________________________ >> > > Libmesh-users mailing list >> > > Libmesh-users@lists.sourceforge.net >> > > https://lists.sourceforge.net/lists/listinfo/libmesh-users >> > >> >> >> -- >> Yuxiang "Shawn" Wang, PhD >> yw...@virginia.edu >> +1 (434) 284-0836 >> >> _______________________________________________ >> Libmesh-users mailing list >> Libmesh-users@lists.sourceforge.net >> https://lists.sourceforge.net/lists/listinfo/libmesh-users >> > -- Yuxiang "Shawn" Wang, PhD yw...@virginia.edu +1 (434) 284-0836 _______________________________________________ Libmesh-users mailing list Libmesh-users@lists.sourceforge.net https://lists.sourceforge.net/lists/listinfo/libmesh-users
