any suggestions to speed up explicit DG updates?

More
4 years 4 months ago #2909 by Guosheng Fu
Dear all,

I have attached a demo code for incompressible two phase flow using hybrid-mixed FEM (RT0) for pressure and DG-P0 for saturation.
Currently, I use the standard DG updates for the saturation conservation law equation.
I believe this default version could be significantly improved, but I don't know what's the best option available now in ngsolve. I would like to get ideas about speeding up this code. Many thanks as always!

Best,
Guosheng
More
4 years 4 months ago #2923 by mneunteufel
Hi Guosheng,

I tried to modify your ForwardEuler by replacing the Matrix-Vector multiplication, but this did not really help.
Code:
def ForwardEuler(sh, dt): rhsS.data = aS.mat * sh #rhsS.data = invm*rhsS fesS.SolveM(rhsS) return sh-dt*rhsS
Also to "inline" it did not bring a lot
Code:
#gfuS.vec.data = ForwardEuler(gfuS.vec, dt) aS.Apply(gfuS.vec,rhsS) gfuS.vec.data -= dt*invm*rhsS

However, the Compile flag made a huge difference! It is important to compile everything, not only a part of it
Code:
if kdg>0: aS += (-uh*fw(ss)*grad(tt)).Compile(true_compile, wait=True)*dx ... aS += (uh*n*Fhat*tt).Compile(true_compile, wait=True)*dx(element_boundary=True)
With this I got a speed up of a factor 3!

It is also possible to speed up your "fast" part even more. Instead of calling mat.Invert every time you can "update" the inverse:
Code:
a.Assemble() inv = a.mat.Inverse(fesP.FreeDofs(True),inverse="sparsecholesky") while time < tend-dt/2 : if stepS%nsubcycles == 0: a.Assemble() inv.Update() ... gfuP.vec.data = inv * rhsP

Then the symbolic factorization is done only once and not every step. This improved the "fast" part from 2.4 to 1.5 seconds on my machine.

Best
Michael
More
4 years 4 months ago #2925 by Guosheng Fu
Hi Michael,

Thanks about the two tricks, they are helpful.
I also get roughly a factor of 3 speedup with the correct use of the Compile flag, and a factor of 2 for the sparse Cholesky inverse update.

Changing element_boundary formulation to skeleton-based formulation further helps to speed-up the explicit DG solver slightly.

I am now thinking about using implicit DG solver for saturation. But Newton solver fails since .Other() only works for proxyfunctions, not gridfunctions. Do you have any suggestions on this issue?

Best,
Guosheng
Time to create page: 0.103 seconds