Forum Message

 

 

We have moved the forum to https://forum.ngsolve.org . This is an archived version of the topics until 05/05/23. All the topics were moved to the new forum and conversations can be continued there. This forum is just kept as legacy to not invalidate old links. If you want to continue a conversation just look for the topic in the new forum.

Notice

The forum is in read only mode.

Question on residual update costs

More
4 years 2 weeks ago #2548 by Guosheng Fu
Question on residual update costs was created by Guosheng Fu
Hello,

I am trying to speed up my HDG Naver Stokes solver and need some help on technical details.
In my 2D lid driven cavity problem, I am using a sparse cholesky factorization to solve the linear system,
and all integration terms are already appended with the command
Code:
.Compile(True,True)
.

I found that for a particular problem (k=2 on 128 X 128 grid), the cost of residual update is about twice more expensive than the linear system solver. So I am looking into alternative ways to speed up the residual updates. The residual term contains the nonlinear convection and linear viscosity operators.
I tried two ways for the viscosity operator updates:
Code:
visc = LinearForm(V) visc += Grad(uh)*Grad(v)*dx ... visc.Assemble()
and
Code:
visc =BinearForm(V, nonassemble=True) visc += Grad(u)*Grad(v)*dx ... res = visc.mat*uh.vec
Surprisingly, the second BilinearForm approach is about 4 times faster than the LinearForm approach. Why is this the case? I thought the two approaches shall be the same (is it faster to take derivatives of proxi functions than grid functions??).

Also, attached is the snippet for my residuals:
Code:
#### diffusion update: very expensive visc = BilinearForm(fes, nonassemble=True) visc += (-nu*InnerProduct(gradu,gradv)).Compile(True,True)*dx visc += (nu*(gradu*n*tang(v-vhat)+gradv*n*tang(u-uhat) -penalty*tang(u-uhat)*tang(v-vhat))).Compile(True,True)*dx(element_boundary=True) conv = BilinearForm(fes, nonassemble=True) conv += (-InnerProduct(Grad(u)*u, v)).Compile(True,True)*dx vminus = IfPos(u*n, v.Other(), v) conv += u2*n*(u-u.Other())*vminus.Compile(True,True)*dx(skeleton=True) .... res.data = visc.mat*gfu.vec + conv.mat*gfu.vec
Are there immediate ways to further speed up this step? Thank you.

Best,
Guosheng
More
4 years 2 weeks ago #2552 by schruste
Replied by schruste on topic Question on residual update costs
Hi Guosheng,

Why don't you assemble the viscosity bilinear form and apply the matrix? If you do so, take the nonsymmetric storage of the matrix (it has a faster matrix-vector multiplication than the storage-saving symmetric storage).

For the convection operator there are different tricks around. If you have an Hdiv space in your HDG discretization it is typically more efficient to project into an L2-space, apply the convection there and project the functional back. I recommend having a look at the two model templates for NavierStokes which you can find here: github.com/ngsolve/modeltemplates

Best,
Christoph
More
4 years 1 week ago #2554 by Guosheng Fu
Replied by Guosheng Fu on topic Question on residual update costs
Thanks, Christoph. This is helpful. I thought the matrix version shall be slower than assembly on the fly... anyway, I ended up not updating the viscous part on the right hand side.

Also, I (re)found a bug in my code when applying Dirichlet BC conditions for static condensed linear system and found out the remedy in bvp.py (just put it here in case I forgot about it later):
we shall first set zero all internal degrees of freedom before applying the boundary conditions which totally make sense.

# zero local dofs (THIS IS IMPORTANT FOR CONDENSATION!!!!)
innerbits = fes.FreeDofs(False) & ~fes.FreeDofs(True)
Projector(innerbits, False).Project(gfu.vec)
rhs.data -= a.mat*gfu.vec # add dirichlet BC

Time to create page: 0.126 seconds

Login