using the symmetric stress tensor in the Stokes part
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
TensorCorrect = -nu*InnerProduct((grad(u).trans)*n,v)
I have added this to my bilinear form using
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,
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:
- Check if the problem assembles if you remove the definedon=… in the integral
- Attach a minimal example
removing the definedon= flag did not solve the problem, and I have attached a minimal example
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,
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.
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 )
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 ?
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.
InnerProduct(grad(u)+grau(u).trans, grad(v)) is symmetric as
InnerProduct(grad(u)+grau(u).trans, grad(v)) = InnerProduct(grad(u)+grau(u).trans, grad(v).trans).