Integration of the interpolated initial condition

Linear Algebra
1

Hello, I want to calculate the linear form L manually, which contains interpolated initial condition. I have faced some difficulties which related to integration, code:

from dolfin import *
import matplotlib.pyplot as plt
import numpy as np
import sympy as sym

mesh = Mesh('untitled.xml')
parameters["reorder_dofs_serial"] = False

plot(mesh)
plt.show()

initial_condition = Expression('5*x[0]', degree=2)						#Constant(1)

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

u_1 = interpolate(initial_condition, V)

print('fenics interpolation:', u_1.vector().get_local())

u = TrialFunction(V)
v = TestFunction(V)
	
L = u_1*v*dx

b = assemble(L)

print('fenics b vector:', b.get_local())


################## integration ####################
n, e = sym.symbols('n e')

#basis functions
phi_local0 = e
phi_local1 = n
phi_local2 = 1-e-n

#element 0, points 0, 4, 7
i = [0, 1]
j = [0.5, 1]
k = [0, 0.5]

detJacobian = abs(np.linalg.det(np.array([[1, i[0], i[1]], [1, j[0], j[1]], [1, k[0], k[1]]])))

b_new = sym.zeros(1, 1)
b_new = sym.integrate(u_1.vector().get_local()[0] * phi_local0  * detJacobian, (n, 0, 1-e), (e, 0, 1))

print('\nb vector (index 0):', b_new)


############################### dofs ###################################
gdim = mesh.geometry().dim()
print('\nDof:', V.dofmap().dofs())

dofmap = V.dofmap()
dofs = dofmap.dofs()
	
dofs_x = V.tabulate_dof_coordinates().reshape((-1, gdim))

x, y = [], []

#dof of point
for dof, dof_x in zip(dofs, dofs_x):
	print('Point', dof, ':', dof_x)
	x.append(dof_x[0]), y.append(dof_x[1])

untitled.xml

<?xml version="1.0" encoding="UTF-8"?>

<dolfin xmlns:dolfin="http://www.fenicsproject.org">
  <mesh celltype="triangle" dim="2">
    <vertices size="9">
      <vertex index="0" x="0.0000000000000000e+00" y="1.0000000000000000e+00" z="0.0000000000000000e+00"/>
      <vertex index="1" x="0.0000000000000000e+00" y="-0.0000000000000000e+00" z="0.0000000000000000e+00"/>
      <vertex index="2" x="1.0000000000000000e+00" y="0.0000000000000000e+00" z="0.0000000000000000e+00"/>
      <vertex index="3" x="1.0000000000000000e+00" y="1.0000000000000000e+00" z="0.0000000000000000e+00"/>
      <vertex index="4" x="4.9999999999869210e-01" y="1.0000000000000000e+00" z="0.0000000000000000e+00"/>
      <vertex index="5" x="1.0000000000000000e+00" y="5.0000000000205946e-01" z="0.0000000000000000e+00"/>
      <vertex index="6" x="5.0000000000205946e-01" y="0.0000000000000000e+00" z="0.0000000000000000e+00"/>
      <vertex index="7" x="0.0000000000000000e+00" y="4.9999999999869210e-01" z="0.0000000000000000e+00"/>
      <vertex index="8" x="5.0000000000037581e-01" y="5.0000000000037570e-01" z="0.0000000000000000e+00"/>
    </vertices>
    <cells size="8">
      <triangle index="0" v0="0" v1="4" v2="7"/>
      <triangle index="1" v0="7" v1="4" v2="8"/>
      <triangle index="2" v0="7" v1="8" v2="1"/>
      <triangle index="3" v0="1" v1="8" v2="6"/>
      <triangle index="4" v0="4" v1="3" v2="8"/>
      <triangle index="5" v0="8" v1="3" v2="5"/>
      <triangle index="6" v0="8" v1="5" v2="6"/>
      <triangle index="7" v0="6" v1="5" v2="2"/>
    </cells>
  </mesh>
</dolfin>

After interpolating initial condition we get vector of this interpolated function on mesh. And in assembling vector b, b_new = sym.integrate(u_1.vector().get_local()[0] * phi_local0 * detJacobian, (n, 0, 1-e), (e, 0, 1)) what we need to integrate an vector, an element of the vector? I can’t figure out. Thanks for your help

2

Forgot the basis function in initial condition term, then, there is no need integrate a vector, it has to be multiplied by mass matrix