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.

Neumann BC correction with static condensation

More
5 years 2 months ago #1416 by hvwahl
Hi

using the symmetric stress tensor in the Stokes part
Code:
Stokes = nu*InnerProduct(grad(u)+grad(u).trans,grad(v)) - p*div(v) - q*div(u)
of my bilinear form for a flow problem, I need to correct the homogeneous Neumann BC with
Code:
TensorCorrect = -nu*InnerProduct((grad(u).trans)*n,v)
I have added this to my bilinear form using
Code:
a += SymbolicBFI(form = TensorCorrect, VOL_or_BND=BND, skeleton=True, definedon=mesh.Boundaries("outflow"))
Unfortunately, if I attempt to also use static condensation, ( BilinearForm(..., condense=True) ), assembly of the matrix failes with the following output:

RuntimeError: SparseMatrixTM::AddElementMatrix: illegal dnumsin Assemble BilinearForm 'biform_from_py'

Adding dgjumps=True to the FESpace allows me to assemble the matrix, however I would like to avoid this loss of sparsity. Is there a better way to enable the outflow correction with static condensation?

Thank you in advance and best wishes,
Henry
More
5 years 2 months ago #1417 by schruste
Hi Henry,

Although I don't see why this correction should make sense (shouldn't the stress be the quantity that you apply your boundary condition on?), there might be two things that you can do for helping me spotting the problem:
1. Check if the problem assembles if you remove the definedon=... in the integral
2. Attach a minimal example :)

Best,
Christoph
More
5 years 2 months ago #1419 by hvwahl
Hi Christoph,

removing the definedon= flag did not solve the problem, and I have attached a minimal example :)

File Attachment:

File Name: MWE_FullTe...dense.py
File Size:2 KB


When the boundary condition is applied to the full symmetric stress tensor, the solution is allowed to "fan out" at the outflow boundary. We on the other hand, want the solution to flow out parallel to the channel walls, in order to imitate an "infinite channel" at the outlet, and therefore only consider the "standard" Neumann do-nothing condition, leaving this boundary integral left over.

Thank you for your help and best wishes,
Henry
More
5 years 2 months ago #1423 by joachim
Hi Henry,

without static condensation and also without dgjumps your example is working.

The loops over elements and facets are different. Both compute element matrices, and assemble them independently. This cannot fit together with static condensation.

workaround:
You can use SymbolcBFI (..., element_boundary=True). These terms are computed within the element loop. define a FacetFESpace - Gridfunction indicator, set it 1 at the outflow boundary, and 0 else, and use it within the integrator:

SymbolicBFI ( form = indicator * TensorCorrect )


BUT:
I am also not convinced by this modification. You give up symmetry of the matrix, and stability needs also an extra argument. Is there a symmetric and elliptic version ?

What about the alternative to set homogeneous Dirichlet for the tangential velocity at the outlet ?

Best, Joachim
More
5 years 2 months ago #1440 by hvwahl
Hi Joachim,

thank you for you help! Your workaround has worked.

Am not aware of a symmetric and elliptic version. The symmetry of the matrix is already lost by using the symmetric tensor InnerProduct(grad(u)+grau(u).trans, grad(v)). However I'll look into, whether it's easier, to keep the Stokes part in the standard formulation, and to use a separate form to compute the full Cauchy stress tensor for drag, torque etc.

Best wishes,
Henry
More
5 years 2 months ago #1441 by schruste
Hi Henry,
Code:
InnerProduct(grad(u)+grau(u).trans, grad(v))
is symmetric as
Code:
InnerProduct(grad(u)+grau(u).trans, grad(v)) = InnerProduct(grad(u)+grau(u).trans, grad(v).trans)
.

Best, Christoph
Time to create page: 0.193 seconds