Hi, I am modeling a 3D magnetostatic problem in which I want to describe the current density in a circular coil using Nedelec elements. The current density is circular and I define it myself. When I visualize the it by interpolating from the original Nedelec to DG space, I see something very irregular. What is more worrisome, even the Jz component is very strongly-nonzero at points when it is originally defined to be identically zero everywhere. Below is a screenshot of the paraview visualization of the magnitude distribution (min, max should be 0,1 ideally): [image] The function I use to define the current density gives the ideal distribution whose magnitude looks like this: [image] I am unsure where the problem lies, if at all. Is this is a bug/inaccuracy in the interpolation process in dolfinx, how things can be best fit into these mathematical functions or Paraview is doing a poor job at visualization? Below is an MWE that generates the Js.xdmf file which can be visualized in Paraview. It should be run in dolfinx 0.6.0 because the latest release somehow does not work with the last part of saving as xdmf file import dolfinx, gmsh, mpi4py, ufl import numpy as np import matplotlib.pyplot as plt from dolfinx.io import gmshio cyl_rad = 0.1, 0.2 cyl_height = 0.25 fov_rad = 0.5 #%% create a cylindrical shell inside a box def create_geom(disp=False): gmsh.initialize() gmsh.clear() gdim = 3 occ = gmsh.model.occ # inner cylinder obj_cyl_disk = occ.addDisk(0, 0, 0, cyl_rad[0], cyl_rad[0]) _, obj_cyl1, _ = occ.extrude([(gdim-1, obj_cyl_disk)], 0, 0, cyl_height) # outer cylinder obj_cyl_disk = occ.addDisk(0, 0, 0, cyl_rad[1], cyl_rad[1]) _, obj_cyl2, _ = occ.extrude([(gdim-1, obj_cyl_disk)], 0, 0, cyl_height) # enclosure obj_fov = occ.addBox(-fov_rad, -fov_rad, 0, 2*fov_rad, 2*fov_rad, fov_rad) final = occ.fragment([(gdim, obj_fov)], [obj_cyl1, obj_cyl2]) occ.synchronize() # tag domains and boundaries all_doms = gmsh.model.getEntities(gdim) for idx, dom in enumerate(all_doms): gmsh.model.addPhysicalGroup(dom[0], (dom[1],), idx+1) all_edges = gmsh.model.getEntities(gdim-1) for idx, edge in enumerate(all_edges): gmsh.model.addPhysicalGroup(edge[0], (edge[1],), edge[1]) gmsh.model.mesh.generate(gdim) gmsh.model.mesh.refine() mesh, ct, ft = gmshio.model_to_mesh(gmsh.model, mpi4py.MPI.COMM_WORLD, 0) if disp: gmsh.fltk.run() gmsh.finalize() return mesh, ct, ft mesh, ct, ft = create_geom(disp=False) # make True for visual inspection # We have a circular coil, define a current density vector inside it def current_density(r): '''define circular current density in the coil which is zero everywhere outside''' R_in, R_out = cyl_rad # read the global variable cyl_rad height = cyl_height # height, read from the global variable x, y, z = r rho = np.sqrt(x**2 + y**2) print("Interpolated coordinate radius min value = {}, max value = {}".format(rho.min(), rho.max())) print("Interpolated coordinates' min, max z value = {}, {}".format(z.min(), z.max())) theta = np.arctan2(y, x) # angle theta, relevant only in the rounded regions Jx = np.full_like(x, 0.0) Jy = np.full_like(x, 0.0) Jz = np.full_like(x, 0.0) # this component will remain zero points_in_coil = (rho >= R_in)*(rho =0)*(z<=height) Jx[points_in_coil] = -np.sin(theta[points_in_coil]) Jy[points_in_coil] = np.cos(theta[points_in_coil]) return Jx, Jy, Jz # for electromagnetic problem, we will need current density on curl elements curl_el = ufl.FiniteElement("N1curl", mesh.ufl_cell(), 2) Omega = dolfinx.fem.FunctionSpace(mesh, curl_el) Js = dolfinx.fem.Function(Omega) Js.interpolate(current_density) # For visualization, interpolate to DG space Omega_dg = dolfinx.fem.VectorFunctionSpace(mesh, ("DG", 2)) Js_dg = dolfinx.fem.Function(Omega) Js_dg.interpolate(Js) # save as xdmf to visualize in paraview with dolfinx.io.XDMFFile(mesh.comm, "Js.xdmf", "w", encoding=dolfinx.io.XDMFFile.Encoding.ASCII) as file: file.write_mesh(mesh) file.write_function(Js_dg) Thanks a lot for comments and suggestions.

You should not use XDMFFile to visualize Nedelec functions. XDMFile only supports (continuous) Lagrange. In your case, with 0.6.0 it interpolates the function into a first order continuous space. You should interpolate the solution into a suitable DG space, and use VTXWriter or VTKFile to store the solution.

Thank you for your feedback! I interpolated to DG but I guess xdmf is still not a suitable way to save. With VTXWriter, what is the correct reader in Paraview? It shows a huge list and so far everything I tried led to errors.

Thank you so much! I also have something to add for Windows users: Paraview gives ADIOS2VTXReader option only in the MPI version. For some reason, the non-MPI version that I was using does not have it. Secondly, if the is an error due to missing msmpi.dll, install MS MPI before installing Paraview’s MPI version.

Inaccurate interpolation of function

General

dokken December 3, 2023, 3:57pm 4

See Write mixed element function to files - #11 by dokken

Topic		Replies	Views
Export solution to XDMF (Nédélec elements) General	3	428	December 17, 2022
Export functionalities of Dolfinx 0.9 and troubleshooting with the interpolate function dolfinx	3	111	July 9, 2025
N1curl interpolation issue for dipole sources General interpolation	10	546	January 30, 2024
Maximum value in vector and paraview is different dolfinx	12	1429	December 9, 2021
Interpolating Functions on higher order elements General	7	1699	November 7, 2023

Inaccurate interpolation of function

Related topics