difference between matrix-action and form-action
Dear firedrakers, Assume I have assembled a UFL form in the HDiv space: form_uu = [a UFL form] mat_uu = assemble(form_uu,bcs=bcs) then shouldn’t the following two experssions fo exactly the same thing? version A: with v.dat.vec as y: with u.dat.vec_ro as x: mat_uu.M.handle.mult(x,y) version B: v = self._mutilde._form_action(u) as fas as I can tell, only version B gives the correct result in my case. Thanks, Eike
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 On 27/04/15 12:04, Eike Mueller wrote:
Dear firedrakers,
Assume I have assembled a UFL form in the HDiv space:
form_uu = [a UFL form] mat_uu = assemble(form_uu,bcs=bcs)
This assembles an operator with bcs applied, useable in a linear solver. The rows /and/ columns corresponding to the BC nodes are the identity.
then shouldn’t the following two experssions fo exactly the same thing?
version A: with v.dat.vec as y: with u.dat.vec_ro as x: mat_uu.M.handle.mult(x,y)
So this computes I_{Gamma} + A_{Omega \ Gamma}, i.e. on the boundary nodes y is x elsewhere it's Ax.
version B: v = self._mutilde._form_action(u)
This computes the action of A on u /without/ boundary conditions. What are you expecting to see?, Ax on the non-boundary condition nodes, bc-value on the boundary values?
as fas as I can tell, only version B gives the correct result in my case.
Lawrence -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 iQEcBAEBAgAGBQJVPhtOAAoJECOc1kQ8PEYv+VAIAMqhJa+dIZSlSu20g2BSecv6 WFTlDq/gyTu1kQs3mYrwrPTzAUbLt5jMjZMZdUdQVyXnQC4H1i3zBJch5z9smMmT SDL6zIn120XYK+CjAD8snbpRLU/DPvdJ9DMo0QQ1c9swUSLdLxFEnxYuB2OQXqKI nDF2c0T13+uiYQe9m7DzEHMEmN0TB4NZvez8X/6xPb/2lp0U1wwJDxU/7+nk9DDw Qe5RIf5VVzUg21RhtvDAdZqp4AylBPq7GEpUwmcetv4PqYnjgcQMOYQaZGUW1AuC qxVMP29VsjIUCjUlyE9KL7oBBshBjwci4icpkpgP8DC7YnuQrGvqWrgLBOIMHfI= =NtfA -----END PGP SIGNATURE-----
Hi Lawrence, ok, so if I apply the BCs to u before the multiplication and to v after the multiplication I would expect the result to be the same in both versions, but it doesn’t look like it. Thanks, Eike
On 27 Apr 2015, at 12:19, Lawrence Mitchell <lawrence.mitchell@imperial.ac.uk> wrote:
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 27/04/15 12:04, Eike Mueller wrote:
Dear firedrakers,
Assume I have assembled a UFL form in the HDiv space:
form_uu = [a UFL form] mat_uu = assemble(form_uu,bcs=bcs)
This assembles an operator with bcs applied, useable in a linear solver. The rows /and/ columns corresponding to the BC nodes are the identity.
then shouldn’t the following two experssions fo exactly the same thing?
version A: with v.dat.vec as y: with u.dat.vec_ro as x: mat_uu.M.handle.mult(x,y)
So this computes I_{Gamma} + A_{Omega \ Gamma}, i.e. on the boundary nodes y is x elsewhere it's Ax.
version B: v = self._mutilde._form_action(u)
This computes the action of A on u /without/ boundary conditions. What are you expecting to see?, Ax on the non-boundary condition nodes, bc-value on the boundary values?
as fas as I can tell, only version B gives the correct result in my case.
Lawrence -----BEGIN PGP SIGNATURE----- Version: GnuPG v1
iQEcBAEBAgAGBQJVPhtOAAoJECOc1kQ8PEYv+VAIAMqhJa+dIZSlSu20g2BSecv6 WFTlDq/gyTu1kQs3mYrwrPTzAUbLt5jMjZMZdUdQVyXnQC4H1i3zBJch5z9smMmT SDL6zIn120XYK+CjAD8snbpRLU/DPvdJ9DMo0QQ1c9swUSLdLxFEnxYuB2OQXqKI nDF2c0T13+uiYQe9m7DzEHMEmN0TB4NZvez8X/6xPb/2lp0U1wwJDxU/7+nk9DDw Qe5RIf5VVzUg21RhtvDAdZqp4AylBPq7GEpUwmcetv4PqYnjgcQMOYQaZGUW1AuC qxVMP29VsjIUCjUlyE9KL7oBBshBjwci4icpkpgP8DC7YnuQrGvqWrgLBOIMHfI= =NtfA -----END PGP SIGNATURE-----
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
Hi, I think I just saw a rounding error, and the results are the same. I hadn’t realised that the numbers I compare are very large, so that an absolute difference of 1.E-5 is rounding error. I will simply normalise the numbers by the maximum before doing the comparison. Thanks, Eike
On 27 Apr 2015, at 12:26, Eike Mueller <E.Mueller@bath.ac.uk> wrote:
Hi Lawrence,
ok, so if I apply the BCs to u before the multiplication and to v after the multiplication I would expect the result to be the same in both versions, but it doesn’t look like it.
Thanks,
Eike
On 27 Apr 2015, at 12:19, Lawrence Mitchell <lawrence.mitchell@imperial.ac.uk> wrote:
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1
On 27/04/15 12:04, Eike Mueller wrote:
Dear firedrakers,
Assume I have assembled a UFL form in the HDiv space:
form_uu = [a UFL form] mat_uu = assemble(form_uu,bcs=bcs)
This assembles an operator with bcs applied, useable in a linear solver. The rows /and/ columns corresponding to the BC nodes are the identity.
then shouldn’t the following two experssions fo exactly the same thing?
version A: with v.dat.vec as y: with u.dat.vec_ro as x: mat_uu.M.handle.mult(x,y)
So this computes I_{Gamma} + A_{Omega \ Gamma}, i.e. on the boundary nodes y is x elsewhere it's Ax.
version B: v = self._mutilde._form_action(u)
This computes the action of A on u /without/ boundary conditions. What are you expecting to see?, Ax on the non-boundary condition nodes, bc-value on the boundary values?
as fas as I can tell, only version B gives the correct result in my case.
Lawrence -----BEGIN PGP SIGNATURE----- Version: GnuPG v1
iQEcBAEBAgAGBQJVPhtOAAoJECOc1kQ8PEYv+VAIAMqhJa+dIZSlSu20g2BSecv6 WFTlDq/gyTu1kQs3mYrwrPTzAUbLt5jMjZMZdUdQVyXnQC4H1i3zBJch5z9smMmT SDL6zIn120XYK+CjAD8snbpRLU/DPvdJ9DMo0QQ1c9swUSLdLxFEnxYuB2OQXqKI nDF2c0T13+uiYQe9m7DzEHMEmN0TB4NZvez8X/6xPb/2lp0U1wwJDxU/7+nk9DDw Qe5RIf5VVzUg21RhtvDAdZqp4AylBPq7GEpUwmcetv4PqYnjgcQMOYQaZGUW1AuC qxVMP29VsjIUCjUlyE9KL7oBBshBjwci4icpkpgP8DC7YnuQrGvqWrgLBOIMHfI= =NtfA -----END PGP SIGNATURE-----
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
-----BEGIN PGP SIGNED MESSAGE----- Hash: SHA1 On 27/04/15 13:14, Eike Mueller wrote:
Hi,
I think I just saw a rounding error, and the results are the same. I hadn’t realised that the numbers I compare are very large, so that an absolute difference of 1.E-5 is rounding error. I will simply normalise the numbers by the maximum before doing the comparison.
Aha, I think this will depend on what your boundary condition /value/ is. Here's a bit more explanation: Consider splitting the function space V = V_0 \osum V_g where V_0 are the basis functions that vanish on the boundary nodes and V_g are the basis functions that vanish on the interior. So when I assemble my operator it splits into: [ A_00, A_0g; A_g0, A_gg ] Now, A_g0 is zero (due to choice of basis functions) and A_gg is the identity, giving: [ A_00, A_0g; 0 , I ] When we assemble the operator for solving a system of equations, we /also/ throw away A_0g giving: [ A_00, 0; 0 , I ] And forward substitute A_0g onto the right hand side, since it acts on the boundary condition nodes, for which we know the correct solution value. However, consider some (known) Function u, we wish to compute the action of the operator on u. If we have the assembled operator (with bcs) then we end up with: Au = A_00 u_0 + I u_g where u_g is the value on the bc nodes, and u_0 on the complement. So to have Au satisfy the boundary conditions we can /either/ apply the bcs before or after, since the application of the bc operator is the identity. If you're computing the action: assemble(action(a, u), bcs=bcs) Then something different happens: Au|bcs = (A_00 u_0 + A_0g u_0 + A_gg u_g)|bcs So we first assemble the action of the operator /without/ bcs and then apply the bc values. Now the answer will be correct on the bc nodes (it will contain the bc values) but will contain an additional A_0g u_0 term on some of the non-bc nodes (that overlap with the bc basis functions). To get the correct answer, we need to kill the A_0g u_0 term, to do this, we can apply homogeneous bcs to u, compute the action, and then apply the correct bc values to the result. Like so: hbc = homogenize(bc) hbc.apply(u) Au = assemble(action(a, u), bcs=bc) Cheers, Lawrence -----BEGIN PGP SIGNATURE----- Version: GnuPG v1 iQEcBAEBAgAGBQJVPjEHAAoJECOc1kQ8PEYvDEUIAL6+30Xer/JQsPJpHf+UjCQS ZbaIdfoRc63Gd7WyQr3MVqqgQpOqyNABEWDXunS7vWsLIV1l+5lNHtkhG0Fe58NO L0skC+tAetOiIRswev1FZgrhvIBcLcvYePqbhjr7bCGGZGEYfGajErwv/87iydx9 fQL1TO27bAAH3CRaQHJ1mjLOYsfw1PeEpmh+lFQtKpCDKyuvNeQfGtV49kGYl869 2ERZJngUBaVlBIqLaO67mS094Un+6XHAeTcD3Ok5D20Q8pKd4vGxUT7wuXRSP2ji ZhDZY40xuAXvpXQBQbGxoOy59lQzxMF0GSLpCxrtB4L+0C97SMf5N8lnPsOlY3o= =Yxo8 -----END PGP SIGNATURE-----
participants (2)
-
Eike Mueller
-
Lawrence Mitchell