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.

Inconsistency between Apply() and Assemble() when using .Other() in BND integral

More
6 years 7 months ago #136 by alexschlueter
Hi,

it seems like there is an inconsistency between the two methods of a BilinearForm in a skeleton boundary integral:
  • Apply(): all occurences of proxy.Other() are set to the same values as the proxy itself
  • Assemble(): all occurences of proxy.Other() are set to zero

Here is an example:
Code:
from ngsolve import * from netgen.geom2d import unit_square mesh = Mesh(unit_square.GenerateMesh()) Draw(mesh) fes = L2(mesh, flags={'dgjumps': True}) v = fes.TrialFunction() w = fes.TestFunction() a = BilinearForm(fes) a += SymbolicBFI(v.Other()*w, BND, skeleton=True) g = GridFunction(fes) g.Set(CoefficientFunction(1)) app = GridFunction(fes) ass = GridFunction(fes) a.Apply(g.vec, app.vec) print(app.vec) a.Assemble() ass.vec.data = a.mat*g.vec print(ass.vec) print(Integrate((app-ass)*(app-ass), mesh))

Output:
Code:
2 2 0 0 3.999999999999984

I believe the relevant lines in symbolicintegrator.cpp are here for Apply() and here for Assemble() .

Is this a bug?

Alex
More
6 years 7 months ago #144 by joachim
Hi Alex,

you can get indeed different results from Apply and Assemble + matrix times vector:

The Appyl computes the nonlinear operator (which includes also an affine linear operator from inhomogeneous Dirichlet data). Here the class-name BilinearForm might be misleading, since it can be non-linear in the trial-argument.

The Assemble can only compute a matrix, and assumes you are really providing a bilinear-form (linear in both arguments). AssembleLinearized computes the linearization at a given point.

What is the u.Other() doing for skeleton-BND integrals ?
It uses boundary values only if you give a bnd=... parameter, otherwise it's using the value from inside (and thus you got the non-zero result).

This modification gave me identical results (i.e. 0):
Code:
a += SymbolicBFI(v.Other(bnd=0)*w, BND, skeleton=True)

Joachim
More
6 years 7 months ago #145 by alexschlueter
Hi Joachim,

The Assemble can only compute a matrix, and assumes you are really providing a bilinear-form (linear in both arguments)

Isn't the expression v.Other()*w bilinear anyway, regardless of whether the inner value or 0 is used for v.Other()? So I would have expected Assemble() and Apply() to do the same thing.

What is the u.Other() doing for skeleton-BND integrals ?
It uses boundary values only if you give a bnd=... parameter, otherwise it's using the value from inside (and thus you got the non-zero result).

Apply() uses the inside value, Assemble() always uses zero.

Using bnd=0 gives identical results, but does that mean Other() without argument is allowed to be inconsistent? Another example would be
Code:
a += SymbolicBFI(v*v.Other(2)*w, BND, skeleton=True)
which is also bilinear (on the boundary), uses the bnd argument, but gives different results.
That said, nobody would actually write v.Other(2) instead of 2, and I don't have a real-world example right now, so maybe this issue is not all that important :)

Alex
More
6 years 7 months ago #146 by joachim

Isn't the expression v.Other()*w bilinear anyway, regardless of whether the inner value or 0 is used for v.Other()? So I would have expected Assemble() and Apply() to do the same thing.


you are right, it's a bilinear-expression, but both cases are artificial: why should one use Other() to get the internal value, or to get the 0-function ? and you are right, Assemble and Apply should be consistent in this artificial case, I think setting it to 0 is the better choice (i.e. change the behaviour of Apply).


Currently, there is indeed a real application for the last example you gave. There is an intended difference between a1 and a2:
Code:
vals = CoefficientFunction( [3,7,15] ) a1 += SymbolicBFI(v*v.Other(bnd=vals)*w, BND, skeleton=True) a2 += SymbolicBFI(v*vals*w, BND, skeleton=True)

a1 uses the MappedIntegrationRule (mir) on the boundary, while a2 uses the mir of the volume. The mir contains also the region index, either material index or boundary condition index. So different elements of the vals - CF are chosen.

As soon as we have the .Other() also for ordinary CoefficientFunctions, we don't need this tricky construct via the proxy anymore.
Time to create page: 0.180 seconds