ValueError: Shape Mismatch when Comparing Interpolated Solution to Numerical

#1

I am solving a set of coupled Maxwell’s equations and got the system working. To calculate the error, I use analytical expressions for the two fields and am trying to compare those against my numerical results at degrees of freedom (d.o.f.).

I am trying to evaluate at the d.o.f. because I am using Nedelec elements of the 1st kind, and so think that I cannot the following function (which works perfectly) does not give the complete picture:

``````compute_vertex_values()
``````

However, it turns out that I get the following error when I try to compute the error. I do not understand why, because I am using function_space from the solution to define the space I want to interpolate in.

``````ValueError: operands could not be broadcast together with shapes (1115,) (2230,)
``````

It appears as if the interpolated data is half that of the numerical solution every time. Could someone help me see why? I have seen a few posts on here that discuss this issue, but none of them seem to work. For example, /qa/9444/interpolation-on-mixed-functionspace-shapes-not-matching/?show=9447#a9447 suggests using FunctionAssigner() instead of split(), but my implementation of the attempt did not change anything (perhaps I did it wrong). Here is another related post: /qa/14036/initial-condition-in-mixed-function-spaces/.

Here is a minimal working example:

``````from dolfin import *
import numpy as np

[f, pie] = (1, np.pi)
[omega, k] = (Constant(2*pie*f), 2)
[eps, mu] = (Constant(1.0), Constant(1.0))

[Nx, Ny, Nz] = (5, 5, 5)
mesh = UnitCubeMesh(Nx, Ny, Nz)
P1 = FiniteElement("N1E", tetrahedron, 1) # FiniteElement("N1curl", mesh.ufl_cell(), 1) nedelec edge elements
TH = P1*P1
V = FunctionSpace(mesh, TH)

# Define test functions
test_E, test_H = TestFunctions(V)

# Define functions for velocity and concentrations
u = Function(V)

# Split system functions to access components
trial_E, trial_H = split(u)

F = omega * inner(eps*trial_E, test_H) * dx + omega * inner(mu * trial_E, test_E) * dx

solve(F == 0, u)#, bcs) # solve(a == L, u, bc). NO bc rn.
_u_1, _u_2 = u.split()

# Analytical solution:
EE = Expression(("cos(2*piee*kk*x[2])", "cos(2*piee*kk*x[2])", "0"), piee = pie, kk = k, degree=2)
HH = Expression(("-cos(2*piee*kk*x[2])", "cos(2*piee*kk*x[2])", "0"), piee = pie, kk = k, degree=2)

# Interpolate analytical result:
V1 = _u_1.function_space().collapse()
V2 = _u_2.function_space().collapse()
exact_E = interpolate(EE, V1)
exact_H = interpolate(HH, V2)

nodal_values_EE = exact_E.vector().get_local()
nodal_values_HH = exact_H.vector().get_local()

# Compute error
error_max_EE = (nodal_values_EE - _u_1.vector().get_local()).max()
error_max_HH = (nodal_values_HH - _u_2.vector().get_local()).max()
print('error_max_EE  =', error_max_EE)
print('error_max_HH  =', error_max_HH)``````
#2

Based off /qa/13269/scalar-mixed-function-space-problem/, the solution was to set `deepcopy=True` via `u.split(True)`.

1 Like