Dear all: I am new to fenics, and I wonder if it’s possible to find the gradient w.r.t the MeshFunction? More specified: For a given mesh, with stored simulation result from external numerical model, we denote the f_{i} \in R be the value of the i-th mesh center. So we have no explicit expression for f but only mesh-center value. Consider for the nonlinear poisson problem (g is known function, boundary setting is same as demo nonlinear poisson case), my question is : how we formulate the || \nabla f||^{2}, or gradient of f in fenics, with only information of the mesh and the cell center value ? \nabla \cdot \left( {||\nabla f|{|^2}\nabla u} \right) = g Tried solution: [1] In demo case ‘tensor weighted poisson’, an example of the meshFunction may solve this problem. However, it reports that '‘MeshFunctionDouble can not be converted to any UFL type’ . Thanks for anyone watched this thread. Best regards, Anthony

An oversight perhaps is using a DG space V^p_\text{DG} of degree p=0. Naturally the broken gradient operator defined on a mesh \mathcal{T}_h = \{\kappa\} is \nabla_h := \nabla|_{\kappa \in \mathcal{T}_h} when applied to a function \nabla_h u_h(x) = 0, \forall u_h \in V^0_\text{DG}. So the trick would be to decide how to interpolate your piecewise constant data to something which varies intra-cell. Or choose some other method for interpolation/projection of your cell-wise data which varies whilst conforming with continuity inter-cell.

Hello dokken, thanks for your kindly reply! As your suggestion, did you mean the code should be like this? f = [val_cell_1, val_cell_2, ..., val_cell_nxn] mesh = UnitSquareMesh(n, n) V = FunctionSpace(mesh, "DG", 0) v = Function(V) x = V.tabulate_dof_coordinates() for i in range(x.shape[0]): # cell_number x dimension v.vector().vec().setValueLocal(i, f[i]) grad_of_v = grad(v) do something Sorry for this kind of reply since something wrong with my docker and i currently not able to verified by myself. Thanks for your help!

You have written sudo-code almost identical to the post I linked you to, so yes, it should do what you expect.

Yeah, the code above almost referenced from your post, just to verfify the code behave as i expect. Thank you for your kindly help! It really helps me a lot. Best regards! Anthony

Hello dokken, I tried previous code but something wired happened, the mesh is uni-square mesh with 64 \times 64 the function represent synthetic mesh value, here f(x,y):=x^{2}+y^{2}, so \nabla{f}=[2x, 2y]^{T} The gradient function defined by following code produces 0 at each point, but it should not import numpy as np from dolfin import * # Toy example n = 64 mesh = UnitSquareMesh(n, n) V = FunctionSpace(mesh, "DG", 0) f = Function(V) x = V.tabulate_dof_coordinates() # Define f := x^2 + y^2 for i in range(x.shape[0]): x_, y_ = x[i] f.vector().vec().setValueLocal(i, np.sum(np.power(x[i], 2))) # Loop to check the gradient and coordinate # Only first 5 res to check Lx, Ly = grad(f) for idx, crds in enumerate(x): if idx < 5: print("x, y, f(x,y)", crds, f(crds)) print("residual at crds", f(crds) - np.sum(np.power(crds, 2))) print("grad at x is", Lx(crds)) print("grad at y is", Ly(crds)) print("\n") else: break and the output x, y, f(x,y) [ 0.01041667 0.00520833] 0.000135633680556 residual at crds 0.0 grad at x is 0.0 grad at y is 0.0 x, y, f(x,y) [ 0.00520833 0.01041667] 0.000135633680556 residual at crds 0.0 grad at x is 0.0 grad at y is 0.0 x, y, f(x,y) [ 0.02604167 0.00520833] 0.000705295138889 residual at crds 0.0 grad at x is 0.0 grad at y is 0.0 x, y, f(x,y) [ 0.02083333 0.01041667] 0.000542534722222 residual at crds 0.0 grad at x is 0.0 grad at y is 0.0 x, y, f(x,y) [ 0.04166667 0.00520833] 0.00176323784722 residual at crds 0.0 grad at x is 0.0 grad at y is 0.0 The function value is correctly assigned to mesh and the function is properly constructed, but obviouly the gradient should not be 0, so what may be the problem? Thanks for your help

That might be the solution,I should check it more precicely and spend some time checking the document of Fenics. Thans for your kindly help! Have a nice day :slight_smile:

[Gradient w.r.t MeshFunction]

General

dokken November 24, 2022, 1:43pm 2

You would have to convert your mesh function into a function in the (“DG”, 0) space, see for instance: How to define different materials / importet 3D geometry (gmsh) - #2 by dokken

Topic		Replies	Views
Gradient of basic functions evaluated in each point of the mesh General mesh	3	266	June 7, 2024
Calculating gradient from mesh General mesh	3	432	April 5, 2024
How to get the first order derivative of the solution and store it in a numpy array I/O dolfinx	23	1158	April 22, 2025
Derivatives, gradients etc. in Dolfinx Uncategorized	1	1301	September 19, 2022
How can I interpolate a UFL nabla_grad/grad to a function dolfinx	3	867	October 15, 2022

[Gradient w.r.t MeshFunction]

Related topics