Couldn't map 'v_0' to a float, returning ufl object without evaluation

General
1

Hello, I was prompted when calculating the following code:

Couldn’t map ‘v_0’ to a float, returning ufl object without evaluation.

from fenics import *
import numpy as np
import scipy.sparse as sp
import matplotlib.pyplot as plt

kB = Constant(8.6e-5)  
r = Constant(0.3)   
Sext = Constant(0.0)    
fai = Constant(1e24)  
lamda = Constant(0.3165e-9)
v0 = Constant(3e13) 
f = 0.02    
n1 = Constant(8.2372e18)
E1 = Constant(0.87) 
n2 = Constant(2.2e18)
E2 = Constant(1.0)
nsolute = Constant(3.802e22)

t_total = 10   
num_steps = 500   
dt = t_total / num_steps    

def D(T):    
    return 4.17e-7*exp(-0.39/(kB*T))

def KrW(T):    
    return 3.2e-15*exp(-1.16/(kB*T))

def KrCuCrZr(T):   
    return 2.9e-14*exp(-1.92/(kB*T))

def S(T):    
    return 12.915609e23*exp(-1.04/(kB*T))

def vm(T):
    return D(T)/(nsolute*lamda**2)

def v1(T):
    return v0*exp(-E1/(kB*T))

def v2(T):
    return v0*exp(-E2/(kB*T))

mesh = Mesh()
with XDMFFile("mesh/mesh.xdmf") as infile:
    infile.read(mesh)

W = FunctionSpace(mesh, 'P', 1)

P1 = FiniteElement('P', triangle, 1)
element = MixedElement([P1, P1, P1])
V = FunctionSpace(mesh, element)

class BoundaryX0(SubDomain):   
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[0] - 0.00) < DOLFIN_EPS

class BoundaryX1(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[0] - 0.01) < DOLFIN_EPS
    
class BoundaryY0(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[1] - 0.00) < DOLFIN_EPS

class BoundaryY1(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[1] - 0.01) < DOLFIN_EPS

bx0 = BoundaryX0()    
bx1 = BoundaryX1()
by0 = BoundaryY0()
by1 = BoundaryY1()

boundary_markers = MeshFunction('size_t', mesh, mesh.topology().dim()-1)   

boundary_markers.set_all(9999)  
bx0.mark(boundary_markers, 0)   
bx1.mark(boundary_markers, 1)
by0.mark(boundary_markers, 2)
by1.mark(boundary_markers, 3)

ds = Measure('ds', domain=mesh, subdomain_data=boundary_markers)
        
T = Function(W)

C = Function(V)
C_n = Function(V) 
c, ct1, ct2 = split(C)
c_n, ct1_n, ct2_n = split(C_n)

v1, v2, v3 = TestFunctions(V)

K1 = KrW(T)*S(T)**2

F = ((c - c_n)/dt)*v1*dx + D(T)*dot(grad(c),grad(v1))*dx - Sext*v1*dx - ((1-r)*fai-KrW(T)*c**2)*v1*ds(3) + (KrCuCrZr(T)*c**2)*v1*ds(2) + (KrW(T)*c**2-K1*f)*v1*ds(0) + (KrW(T)*c**2-K1*f)*v1*ds(1) \
+ ((ct1 - ct1_n)/dt)*v2*dx - vm(T)*c*(n1 - ct1)*v2*dx + v1(T)*ct1*v2*dx \
+ ((ct2 - ct2_n)/dt)*v3*dx - vm(T)*c*(n2 - ct2)*v3*dx + v2(T)*ct2*v3*dx

timeseries_T = TimeSeries('heat_transfer_11s/temperature_series')

t = 0

plt.figure(figsize=(12,8))
for n in range(num_steps):
   
    t += dt
    
    timeseries_T.retrieve(T.vector(), t)
        
    solve(F == 0, C)
    
    plot(C)
    
    C_n.assign(C)
2

You need to change v1(T)*ct1*v2*dx to v1*ct1*v2*dx and v2(T)*ct2*v3*dx to v2*ct2*v3*dx in the expression of F. So the modified F is -

F = ((c - c_n)/dt)*v1*dx + D(T)*dot(grad(c),grad(v1))*dx - Sext*v1*dx \
- ((1-r)*fai-KrW(T)*c**2)*v1*ds(3) + (KrCuCrZr(T)*c**2)*v1*ds(2) \
+ (KrW(T)*c**2-K1*f)*v1*ds(0) + (KrW(T)*c**2-K1*f)*v1*ds(1) \
+ ((ct1 - ct1_n)/dt)*v2*dx - vm(T)*c*(n1 - ct1)*v2*dx + v1*ct1*v2*dx \
+ ((ct2 - ct2_n)/dt)*v3*dx - vm(T)*c*(n2 - ct2)*v3*dx + v2*ct2*v3*dx
3

Thank you for your answer, but I have the following questions:
(1) I have modified v(T) to v according to your reply. I would like to know why I need to do so？
(2) After modifying v(T) to v, the following error is displayed:

too many values to unpack (expected 2)

4

v1, v2 etc are functions in your definition:

They take in an input, and return an expression, thus you need to get the expression (as otherwise the input to v1,v2 is unknown).

I cant find FF at all in the post above.

Please post your updated code, as it is unclear what you have changed here.

5

This is my updated code:

from fenics import *
import numpy as np
import scipy.sparse as sp
import matplotlib.pyplot as plt

kB = Constant(8.6e-5)  
r = Constant(0.3)   
Sext = Constant(0.0)    
fai = Constant(1e24)  
lamda = Constant(0.3165e-9)
v0 = Constant(3e13) 
f = 0.02    
n1 = Constant(8.2372e18)
E1 = Constant(0.87) 
n2 = Constant(2.2e18)
E2 = Constant(1.0)
nsolute = Constant(3.802e22)

t_total = 10   
num_steps = 500   
dt = t_total / num_steps    

def D(T):    
    return 4.17e-7*exp(-0.39/(kB*T))

def KrW(T):    
    return 3.2e-15*exp(-1.16/(kB*T))

def KrCuCrZr(T):   
    return 2.9e-14*exp(-1.92/(kB*T))

def S(T):    
    return 12.915609e23*exp(-1.04/(kB*T))

def vm(T):
    return D(T)/(nsolute*lamda**2)

def v1(T):
    return v0*exp(-E1/(kB*T))

def v2(T):
    return v0*exp(-E2/(kB*T))

mesh = Mesh()
with XDMFFile("mesh/mesh.xdmf") as infile:
    infile.read(mesh)

W = FunctionSpace(mesh, 'P', 1)

P1 = FiniteElement('P', triangle, 1)
element = MixedElement([P1, P1, P1])
V = FunctionSpace(mesh, element)

class BoundaryX0(SubDomain):   
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[0] - 0.00) < DOLFIN_EPS

class BoundaryX1(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[0] - 0.01) < DOLFIN_EPS
    
class BoundaryY0(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[1] - 0.00) < DOLFIN_EPS

class BoundaryY1(SubDomain):
    def inside(self, x, on_boundary):
        return on_boundary and abs(x[1] - 0.01) < DOLFIN_EPS

bx0 = BoundaryX0()    
bx1 = BoundaryX1()
by0 = BoundaryY0()
by1 = BoundaryY1()

boundary_markers = MeshFunction('size_t', mesh, mesh.topology().dim()-1)   

boundary_markers.set_all(9999)  
bx0.mark(boundary_markers, 0)   
bx1.mark(boundary_markers, 1)
by0.mark(boundary_markers, 2)
by1.mark(boundary_markers, 3)

ds = Measure('ds', domain=mesh, subdomain_data=boundary_markers)
        
T = Function(W)

C = Function(V)
C_n = Function(V) 
c, ct1, ct2 = split(C)
c_n, ct1_n, ct2_n = split(C_n)

v1, v2, v3 = TestFunctions(V)

K1 = KrW(T)*S(T)**2

F = ((c - c_n)/dt)*v1*dx + D(T)*dot(grad(c),grad(v1))*dx - Sext*v1*dx - ((1-r)*fai-KrW(T)*c**2)*v1*ds(3) + (KrCuCrZr(T)*c**2)*v1*ds(2) + (KrW(T)*c**2-K1*f)*v1*ds(0) + (KrW(T)*c**2-K1*f)*v1*ds(1) \
+ ((ct1 - ct1_n)/dt)*v2*dx - vm(T)*c*(n1 - ct1)*v2*dx + v1*ct1*v2*dx \
+ ((ct2 - ct2_n)/dt)*v3*dx - vm(T)*c*(n2 - ct2)*v3*dx + v2*ct2*v3*dx

timeseries_T = TimeSeries('heat_transfer_11s/temperature_series')

t = 0

plt.figure(figsize=(12,8))
for n in range(num_steps):
   
    t += dt
    
    timeseries_T.retrieve(T.vector(), t)
        
    solve(F == 0, C)
    
    plot(C)
    
    C_n.assign(C)

[/quote]

6

Sorry, I have already solved the problem, the function v1 and v2 should be in the sign conflict with the test function.

7

Hi dokken, I have another question, how do I add the calculated values c, ct1, ct2 and save them to an XDMFFile file?

8

I assume you mean:

Please just try to write them to file as normal, i.e. XDMFFile("function.xdmf").write(c)

9

What if I want to save the values of c+ct1+ct2 to a file?

10

As they are scalar functions, you can use:

import dolfin

mesh = dolfin.UnitSquareMesh(10, 10)
W = dolfin.FunctionSpace(mesh, 'P', 1)

P1 = dolfin.FiniteElement('P', dolfin.triangle, 1)
element = dolfin.MixedElement([P1, P1, P1])
V = dolfin.FunctionSpace(mesh, element)
u = dolfin.Function(V)

# Manipulate data here
u.interpolate(dolfin.Expression(("x[0]", "x[1]", "3*x[1]"), degree=1))


# For outputting
c_o, c1_o, c2_o = u.split(deepcopy=True)
c_o.vector()[:] += c1_o.vector()[:] + c2_o.vector()[:]
with dolfin.XDMFFile("c_sum.xdmf") as xdmf:
    xdmf.write(c_o)

Please note that this topic is diverging from the original topic.

This topic is also covered in many other posts, so please search the forum before adding more questions to this thread.

11

Thank you, I solved my problem, thank you for your help all the time!