# Using Conditional function inside a Project function

**URL:** <https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967>\
**Category:** Linear Algebra\
**Created:** [August 20, 2020, 6:40am UTC](https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967 "2020-08-20T06:40:34Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![CoE2845](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/coe2845/32/1210_2.png) [@CoE2845](https://fenicsproject.discourse.group/u/CoE2845)\
**Post date:** [August 20, 2020, 6:40am UTC](https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967/1 "2020-08-20T06:40:34Z")

</div>

Hi everyone,

I am currently trying to solve a heat equation problem where I want the values of some constants within my heat equation to change when the temperature at a node of my mesh exceed some value (shown as Tmelting in my code). Currently, I am attempting to implement this via a conditional function:

```auto
assign(kappa, project(conditional(gt(T, Tmelting), kappa_l, kappa_s), Space, solver_type = "mumps", form_compiler_parameters={"cpp_optimize": True, "representation":"uflacs", "quadrature_degree": 2} ))
# kappa_l = value of kappa if T > Tmelt
# kappa_s = value of kappa if T < Tmelt

```

Where kappa, kappa\_l and kappa\_s are all tensors. My aim with this code is to have the value of kappa change from kappa\_s to kappa\_l when T \> Tmelt and have this new value fed back into the weak form for the next time step in my problem. However I am currently getting an error code:

```auto
Traceback (most recent call last):
  File "target3D_surface2.py", line 193, in <module>
    assign(kappa, project(conditional(gt(T, Tmelting), kappa_l, kappa_s), Space, solver_type = "mumps", form_compiler_parameters={"cpp_optimize": True, "representation":"uflacs", "quadrature_degree": 2} ))
AttributeError: 'float' object has no attribute '_cpp_object'

```

Orginially I assumed this was from the floats I had previously used for kappa\_l and kappa\_s instead of the tensors required by the project function, but now I am now sure. Could anybody be able to help with determining the root of the error?

The relevant code is shown below:

```auto
from __future__ import print_function
from fenics import *
import numpy as np
import matplotlib.pyplot as plt
from dolfin import *
#from mshr import *
from ufl import as_tensor
from ufl import Index
import math

parameters["allow_extrapolation"] = True
parameters["form_compiler"]["cpp_optimize"] = True
set_log_level(20)

#Configure for dimensions of system
#----------THERMAL--PROPERTIES--(Solid)--------
kappa_s = 138.0 #conductivity [W/m K] 
c_s = 277.17 #Specific Heat Capacity [J/kg K]
rho_s = 10000.0 #Density [kg/m^3]
tau_T = 0.0 #Temperature Gradient Lag
tau_q = 0.0 #Heat Flux Lag
alpha_s = 5.1162e+7 #absorption coefficient (m^-1)
epsilon_s = 0.05 # material emissivity
#---------------------------------------

#----------THERMAL--PROPERTIES--(liquid)--------
Tmelting = 2896.0
kappa_l = 80.0 #conductivity [W/m K] 
c_l = 350.0 #Specific Heat Capacity [J/kg K]
rho_l = 9330.0 #Density [kg/m^3]
alpha_l = 5.1162e+7 #absorption coefficient (m^-1)
epsilon_l = 0.2 # material emissivity
#---------------------------------------

pi = 3.141592653589793
T_am = 300 #ambient vacuum temperature (K)
T_a4 = T_am**4
sigma = 5.67E-8 # W/(m **2.K** 4) S-B Constant
Pmax = 200 #peak power of laser in W
w = 0.56e-3 #m
Iden = Pmax / (pi*w*w) #Intensity 

Laser = Expression('B*A*Iden*exp(-pow((x[0] - 0), 2)/(w*w)-pow((x[1]-0), 2)/(w*w))',degree=2, A=1, Iden=Iden, w= w, B = 1) #power (w2 localises the z-coordinates)
#----------------------------------------------
data = 'mesh_2mm'
mesh = Mesh()
hdf = HDF5File(mesh.mpi_comm(), data+'.h5', 'r')
hdf.read(mesh, '/mesh', False)
cells = MeshFunction('size_t', mesh, 3)
hdf.read(cells, '/cells')
facets = MeshFunction('size_t', mesh, 2)
hdf.read(facets, '/facets')
#------------------------------------------------

Space = FunctionSpace(mesh, 'P', 1) #define finite element function space, defined via basis functions
VectorSpace = VectorFunctionSpace(mesh, 'P', 1)
da = Measure('ds', domain=mesh, subdomain_data = facets, metadata = {'quadrature_degree': 2}) #area element

dv = Measure('dx', domain=mesh, subdomain_data = cells, metadata = {'quadrature_degree': 2}) #volume element

T = Function(Space)
T0 = Function(Space)
T_init = Expression('Tambient', degree=1, Tambient=300.)
T = project(T_init, Space)
assign(T0,T)

V = TestFunction(Space) # Test Function used for FEA
dT = TrialFunction(Space) # Temperature Derivative
q0 = Function(VectorSpace) # heat flux at previous time step
i = Index()
G = as_tensor(T.dx(i), (i)) #gradient of T
G0 = as_tensor(T0.dx(i), (i)) # gradient of T at previous time step 

#assign initial values for weak form
kappa = kappa_s #conductivity [W/m K] 
c = c_s #Specific Heat Capacity [J/kg K]
rho = rho_s #Density [kg/m^3]
epsilon = epsilon_s # material emissivity
A = epsilon
es = epsilon*sigma

q = as_tensor(dtime/(dtime + tau_q) * (tau_q/dtime*q0[i] - kappa*(1+tau_T/dtime)*G[i] + kappa*tau_T/dtime*G0[i]),(i)) #heat
F = (rho*c/dtime*(T-T0)*V - q[i]*V.dx(i)) * dv + es*(T**4 - T_a4)*V*da - A*Laser*V*da(2) #final form to solve
Gain = derivative(F, T, dT) # Gain will be used as the Jacobian required to determine the evolution of a dynamic system 

problem = NonlinearVariationalProblem(F, T, [], Gain)
solver = NonlinearVariationalSolver(problem)

solver.parameters["newton_solver"]["relative_tolerance"] = 1E-3
solver.parameters["newton_solver"]["absolute_tolerance"] = 1E-4
solver.parameters["newton_solver"]["convergence_criterion"] = "residual"
solver.parameters["newton_solver"]["error_on_nonconvergence"] = True
solver.parameters["newton_solver"]["linear_solver"] = "gmres"
solver.parameters["newton_solver"]["lu_solver"]["symmetric"] = False 
solver.parameters["newton_solver"]["maximum_iterations"] = 100
solver.parameters["newton_solver"]["preconditioner"] = "amg"
solver.parameters["newton_solver"]["relaxation_parameter"] = 1.05
solver.parameters["newton_solver"]["report"] = False

solver.parameters["newton_solver"]["krylov_solver"]["nonzero_initial_guess"] = True 
solver.parameters["newton_solver"]["krylov_solver"]["relative_tolerance"] = 1E-5
solver.parameters["newton_solver"]["krylov_solver"]["absolute_tolerance"] = 1E-6
solver.parameters["newton_solver"]["krylov_solver"]["monitor_convergence"] = False
solver.parameters["newton_solver"]["krylov_solver"]["maximum_iterations"] = 100000

file_T = File('/test/solution.pvd')
t = 0
solver.solve()
kappa = as_tensor(kappa,())
kappa_s = as_tensor(kappa_s,())
kappa_l = as_tensor(kappa_l,())
c = as_tensor(c,())
rho = as_tensor(rho,())
epsilon = as_tensor(epsilon,())
assign(kappa, project(conditional(gt(T, Tmelting), kappa_l, kappa_s), Space, solver_type = "mumps", form_compiler_parameters={"cpp_optimize": True, "representation":"uflacs", "quadrature_degree": 2} ))
file_T << (T,t)
assign(q0, project(q, VectorSpace, solver_type= "bicgstab", form_compiler_parameters={"cpp_optimize": True, "representation": "uflacs", "quadrature_degree": 2} ) )
assign(T0, T)

```

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [August 20, 2020, 7:15am UTC](https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967/2 "2020-08-20T07:15:12Z")

</div>

Please add all definitions of kappa and function spaces, such that one can reproduce the error. Please make it a minimal example, removing the time loop and only include those variables needed to produce the error.

This will make it easier for other people to help you, and is of greater help for other people struggling with the same problem.

---

<div class="post-metadata">

**Author:** ![CoE2845](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/coe2845/32/1210_2.png) [@CoE2845](https://fenicsproject.discourse.group/u/CoE2845)\
**Post date:** [August 20, 2020, 7:41am UTC](https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967/3 "2020-08-20T07:41:57Z")

</div>

Ok, this has been altered accordingly.

---

<div class="post-metadata">

**Author:** ![dokken](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/dokken/32/1560_2.png) [@dokken](https://fenicsproject.discourse.group/u/dokken)\
**Post date:** [August 20, 2020, 8:15am UTC](https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967/4 "2020-08-20T08:15:25Z")

</div>

The problem is that you try to assign a function (the result of the projection), to a float.  
The initial kappa has to be a function.  
This is illustrated in the minimal example below. Please note the following:

1. Your proposed example is not minimal, and contains alot of variables not required to reproduce the error. The example below should be used as a guideline for further inquires.
2. To represent kappa, I am using a DG space, as kappa is discontinuous. If you decide to use a CG space, you will get transfer regions at the switch of the material.

```auto
from dolfin import *
from ufl import as_tensor

kappa_s = 138.0 #conductivity [W/m K]
kappa_l = 80.0 #conductivity [W/m K]

mesh = UnitSquareMesh(10,10)

Space = FunctionSpace(mesh, 'P', 1)
T = project(Expression("x[0]+0.1*x[1]", degree=2), Space)

kappa_space = FunctionSpace(mesh, "DG", 0)
kappa = project(kappa_s,kappa_space)
kappa_s = as_tensor(kappa_s,())
kappa_l = as_tensor(kappa_l,())
Tmelting = 0.5
print(kappa.vector().get_local())
kappa.assign(project(conditional(gt(T, Tmelting), kappa_l, kappa_s), kappa_space))
print(kappa.vector().get_local())

```

---

<div class="post-metadata">

**Author:** ![CoE2845](https://yyz2.discourse-cdn.com/free1/user_avatar/fenicsproject.discourse.group/coe2845/32/1210_2.png) [@CoE2845](https://fenicsproject.discourse.group/u/CoE2845)\
**Post date:** [August 20, 2020, 8:28am UTC](https://fenicsproject.discourse.group/t/using-conditional-function-inside-a-project-function/3967/5 "2020-08-20T08:28:16Z")

</div>

Ah I see. Thank you very much! I will attempt this now
