Question on plotting gradient of a function

Errors
, ,
1

I am trying to plot the gradient of a function I solved to observe the discontinuities in the shape functions. I am getting an error when I run the code.

image
image1423×607 22.6 KB

I made the gradient as dTh and I am trying to plot that. My code is as follows:

from fenics import *

import matplotlib.pyplot as plt

import numpy as np

Generate mesh

domain = IntervalMesh(1000, 0, 1)

SC = FunctionSpace(domain, ‘CG’, 1)

Define boundary and apply BC

def L_bc(x, on_boundary):

return near(x[0], 0) and on_boundary

def R_bc(x, on_boundary):

return near(x[0], 1) and on_boundary

bc = [DirichletBC(SC, Constant(0.), L_bc), DirichletBC(SC, Constant(0.), R_bc)]

Define

#f = Expression(‘exp(x[0]*x[0])’, degree=2)

#C = Expression(‘pow(x[0],3) + 10’, degree=3)

f = Expression(‘exp(-x[0])’, degree=2)

#C = Expression(‘2+x[0]exp(x[0])(x[0]-1)’, degree=2)

C = Expression(‘2*(6*x[0]x[0] - 6x[0] + 1) + x[0]x[0]exp(-x[0])(x[0]-1)(x[0]-1)’, degree=2)

Define test and trial functions

u = TrialFunction(SC)

v = TestFunction(SC)

Weak form

#a = dot(grad(w), grad(T))dx + CwTdx

sol = dot(grad(u), grad(v))dx + vfudx

#L = fwdx

L = vCdx

Compute solution

Th = Function(SC, name=“displacement_linear_1000”)

solve(sol == L, Th, bc)

sigma = grad(u)

x_plt = np.linspace(0, 1, 1000)

Th_plt = np.fromiter(map(Th, x_plt), dtype=np.double)

plt.plot(x_plt, Th_plt)

plt.show()

########################################################

V = VectorFunctionSpace(domain, “DG”, 1)

dTh = project(grad(Th), V)

dTh_plt = np.fromiter(map(dTh, x_plt), dtype=np.double)

plt.plot(x_plt, dTh_plt)

file = File(“displacement_linear_1000.pvd”)

file << dTh

2

Try using XDMFFile and the write_checkpoint function.

If that doesn’t work, I would try

V = FunctionSpace(domain, “DG”, 1)

dTh = project(Th.dx(0), V)

and then use write_checkpoint, as you are in 1D, and then the gradient is the same as d/dx

3

Thanks for the reply. I added the lines you mentioned and am getting the attached image. I ran the simulation for 10, 100 and 1000 nodes. I am expecting to see discontinuities. Are the step functions indicative of the discontinuities?

Results
Results1309×747 51.5 KB

4

Well, yes, as step functions are discontinuous.

As you are using a “CG”-1 function for the primary variable, its derivative is s cell-wise constant.

If you change the degree of the primary variable to 2 you should see slope discontinuities.

5

I understand now, thank you very much for the help.