Hello people,
To make it simple, i have defined a function D_el(x). When i evaluate it with x=C a scalar, i have the expected result. When i evaluate it with x being a field (with x a constant field x=C), i have another (then, wrong) result.
Below is a minimum code that illustrates my problem:
# Impose / is the decimal division (must be at the start of the code)
from __future__ import division
# Dolfin library
from dolfin import *
# Numpy library
import numpy as np
# Current version of FEniCS
print "The current version of FEniCS is",dolfin.dolfin_version()
# Mesh
mesh = BoxMesh(Point(0.0, 0.0, 0.0), Point(10, 10, 10), 10, 10,10)
# Define function
def D_el(Ce,Temp):
C=Ce*1e-6
D=0.00584*exp(-2870/Temp)*((1000*C)**2) - 0.0339*exp(-2920/Temp)*1000*C + 0.129*exp(-3200/Temp)
return(D)
# Test function with scalar
# Input
C=10000
Temp=298.15
# Output
D=D_el(C,Temp)
print D # Correct
# Set function space
V = FunctionSpace(mesh, 'Lagrange', 1)
# Project C on the function space
c_on_V = project(C,V)
plot(c_on_V, title='C',interactive=True) # Value display is correct
# Test function with field
D=D_el(c_on_V,Temp)
plot(D, title='D_on_V',interactive=True) # WRONG !
# Let's get the min and max of the function calculated on the field
D_on_V = project(D,V)
array_=D_on_V.vector().array()
max_=np.max(array_)
print max_ # WRONG !
min_=np.min(array_)
print min_ # WRONG !
I obtain this (all the the three number should be identical)
The current version of FEniCS is 2016.1.0
2.24353357359e-05
Object cannot be plotted directly, projecting to piecewise linears.
-1.6102060254e-05
-1.6102060254e-05
Where is my error?