Replacing FacetNormal produces conj
By replacing the second argument of an inner product with a FacetNormal, one obtains the conjugate of the integral. This is not correct.
Minimal example:
from dolfin import *
from ufl import replace
mesh = UnitSquareMesh.create(10,10, CellType.Type.triangle)
u = Constant((2,0))
v = Constant((1,2))
n = FacetNormal(mesh)
L = inner(u,v)*ds(domain=mesh)
Lun = replace(L, {u: n})
Lvn = replace(L, {v: n})
print(Lun)
print(Lvn)
Produces the following output:
{ (n) : (f_7) } * ds(<Mesh #3>[everywhere], {})
{ conj(((n) : (f_6))) } * ds(<Mesh #3>[everywhere], {})
Comments (10)
-
-
reporter It appears that the conjugate is dependent of the order of declaration of variables. That is fucked up. Minimal example:
from dolfin import * a = Constant((1.0,2.0)) c = Constant((4.0,5.0)) b = Constant((2.0,3.0)) print(inner(a,b)) print(inner(b,a)) print("...") print(inner(b,c)) print(inner(c,b))
Output:
(f_0) : (f_2) conj(((f_0) : (f_2))) ... conj(((f_1) : (f_2))) (f_1) : (f_2)
-
Anyway, there should never be a conjugate because it is built in in "inner", "dot" or ":" and "." respectively already....
-
So
inner(a, b)
gets transformed, later in the day, intosum(a[i]*conj(b[i]), i)
. This is where the conjugate is inserted.To allow for common subexpression elimination. Where you might have
inner(a, b)
andinner(b, a)
, in non-complex UFL these are canonicalised such that the ordering of operands is consistent.Check out what happens on master:
In [2]: a = Constant((1.0,2.0)) ...: c = Constant((4.0,5.0)) ...: b = Constant((2.0,3.0)) ...: ...: print(inner(a,b)) ...: print(inner(b,a)) ...: print("...") ...: print(inner(b,c)) ...: print(inner(c,b)) ...: ...: (w_0) : (w_2) (w_0) : (w_2) ... (w_1) : (w_2) (w_1) : (w_2)
Notice how
inner(a, b)
andinner(b, a)
have the same string representation. This is because they are the same.To allow for this to still occur in complex mode, we need to sort the operands. But, this might result in an invalid
inner
. So, we use the identityinner(a, b) == conj(inner(b, a))
.And now we can sort the operands. This sorting is based on the
count
of a Coefficient, which of course, varies according to the order in which you defined the coefficients.In summary, not a bug, think harder.
-
reporter @wence How do you get this string representation. Using docker with latest dev, i do not get this.
In [4]: from dolfin import * ...: ...: a = Constant((1.0,2.0)) ...: c = Constant((4.0,5.0)) ...: b = Constant((2.0,3.0)) ...: ...: print(inner(a,b)) ...: print(inner(b,a)) ...: print("...") ...: print(inner(b,c)) ...: print(inner(c,b)) ...: ...: (f_9) : (f_11) conj(((f_9) : (f_11))) ... conj(((f_10) : (f_11))) (f_10) : (f_11)
-
I did think hard about this, and I can confirm inner(a,b) and inner(b,a) should have the same string representation, but they don't... that's the problem....
How does your example above resolve constants to "w_", but in my fenics, they are resolved to "f_"?
-
inner(a, b)
andinner(b, a)
cannot have the same representation. And this can't be stressed more. As Lawrence said inner is roughly equivalent tosum(a[i]*conj(b[i]), i)
. This is a decision which has been taken for definition ofinner
. I believe that time to debate about definition ofinner
is over. Please follow the#complex
channel on Slack: https://fenicsproject.slack.com/messages/C93RPKTNC/ -
- changed status to invalid
-
Why should
inner(a, b)
andinner(b, a)
have the same string representation, surely, they should map toa : b
andb : a
(i.e., the argument order is changed).w_
is because I happened to run in firedrake (which labels constant differently), but the UFL issues is the same, I think.I should note that the UFL I was using is a few weeks out of date (so complex is not merged). Hence no
conj
anywhere. -
- changed status to closed
Not a bug, but see https://bitbucket.org/fenics-project/ufl/issues/107 for plausible changes to UFL that would perhaps make things easier for Jørgen and Stephan.
- Log in to comment
Also, the behavior of inner and conjugate is completely random and unclear. Have a look at this:
The output is:
So right now, I see no reason, why apparently "a" should be the first argument and "b" is the second to avoid a conjugation. In any case, the conjugate should be part of "inner" or the ":", so the form/string representation in either case should not get a conjugate, as it is part of the "inner" or the ":" anyway...