On_boundary Clarification

mesh
21

Thanks for answer. But I get this ERROR.

line 55, in <module>
    kappa = K(subdomains, k_0, k_1, degree=0)
********
ufl.log.UFLException: Order "0" invalid for "Lagrange" finite element.

That’s why I write 1 and not 0.

22

Instead of sending in degree, you could send in element=ufl.FiniteElement("DG", mesh.ufl_cell(), 0)

23

Don’t get it clearly. I have to add this stroke and use this element in Function Space:

V = FunctionSpace(mesh, element)

or in stroke with interpolate?

kappa_function = interpolate(kappa, element)

During my attempt, I use element in kappa_function, but I get ERROR:

File "/usr/lib/python2.7/dist-packages/ufl/cell.py", line 359, in as_cell
    error("Invalid cell %s." % cell)
  File "/usr/lib/python2.7/dist-packages/ufl/log.py", line 151, in error
    raise self._exception_type(self._format_raw(*message))
ufl.log.UFLException: Invalid cell [[671 706 716 732]
 [168 197 706 732]
 [679 680 691 742]
 ...,
 [119 491 578 754]
 [119 490 491 754]
 [490 491 610 754]].

Shows that some cell invalid and it makes me confused, why only several cells and not all ('cause in this case I at least start thinking that all my mesh wrong)?

24

Without actually supplying a minimal reproducible example, it is really hard for me to give you any further guidance.

I also note that you are using python 2.7, which was sunset in 2.7.
I would strongly suggest you to use a newer version of Python, see: Sunsetting Python 2 | Python.org

25

This is link on Google Disk where I put my XML file and MED file from Salome.
And this is the whole code

from fenics import *
import numpy as np
import os

k_0 = 1.1
k_1 = 2.1
# Create or using mesh
mesh = Mesh("stripes.xml")

# Define subdomains
tol = 1e-14
class Omega_0(SubDomain):
    def inside(self, x, on_boundary):
        return (0 <= x[0] <= 40 + tol) or (60 - tol <= x[0] <= 100 + tol)

class Omega_1(SubDomain):
    def inside(self, x, on_boundary):
        return (40 - tol <= x[0] <= 60 + tol) or (100 - tol <= x[0] <= 120 + tol)        

# Define conductivity for different subdomains

subdomains = CellFunction('size_t', mesh)
subdomain_0 = Omega_0()
subdomain_1 = Omega_1()
subdomain_0.mark(subdomains, 0)
subdomain_1.mark(subdomains, 1)

# Visualize the subdomains (for illustration purposes)
plot(subdomains, title='Subdomains')

file_s = File('subdomains.xml.gz') << subdomains

# use the values of the mesh function materials to define the variable coefficient kappa 
class K(Expression):
    def __init__(self, subdomains, k_0, k_1, **kwargs):
        self.subdomains = subdomains
        self.k_0 = k_0
        self.k_1 = k_1
    def eval_cell(self, values, x, cell):
        if self.subdomains[cell.index] == 0:
            values[0] = self.k_0
        else:
            values[0] = self.k_1

kappa = K(subdomains, k_0, k_1, degree=1) # degree 1 or 0? 

V = FunctionSpace(mesh, 'DG', 0)

kappa_function = interpolate(kappa, V)

# Save file for checking how we get kappa
file_k = File('kappa_stripes.xml.gz') # 'kappa_stripes.pvd'
file_k << kappa_function

# Get the function values at mesh vertices
kappa_values = kappa_function.compute_vertex_values(mesh)

# Sava values to a txt file
np.savetxt('kappa_values.txt', kappa_values)

Unfortunately, for my simulation I need Python 2.7 and can’t use the new version. In this case, perhaps you can’t help me with it, but could you at least recommend me some literature/link where I can read more about boundary, interpolate in FEniCS and how I can make my own conditions. Thanks in advance.