A problem with the weak form


121
views
0
9 weeks ago by
Nima  

Hi guys,

So I have this chemotaxis term in my integral given by
$\nabla.\left(M\nabla P\right)$∇.(MP)
with the Neumann boundary condition:
$\frac{\partial P}{\partial n}=0$Pn =0
I tried the weak form
$\int\left(\nabla M.\nabla P\right)v\text{ dx}-\int\nabla P.\nabla\left(Mv\right)\text{ dx}$(M.∇P)v dxP.∇(Mv) dx

with the formula
dot(grad(M),grad(P))*v*dx-dot(grad(P),grad(M*v))*dx​
but I receive the error "unsupported operand type(s) for *: 'Measure' and 'Indexed'
Can someone tell me what is wrong?

Here is the complete code:

from dolfin import *
import numpy as np
Domain=Mesh('general_boundary_wenrui+5.xml')
Bulk = MeshFunction('size_t' , Domain , 'general_boundary_wenrui+5_physical_region.xml' ) Boundary = MeshFunction('size_t', Domain , 'general_boundary_wenrui+5_facet_region.xml')

ds = Measure('ds', subdomain_data=Boundary)
dx = Measure('dx',subdomain_data= Bulk)

Dh = 3.93
Dl = 29.89
Dr = 0.205
Dlox = 29.89
Dm = 8.64e-7
Dp = 17.28
kl = 2.35e-4
kh = 5.29e-6
L0 = 9
r0 = 0.26
H0 = 7
alphaH = 1
alphaL = 1
alphaM = 0.2
Lambda = 10
LambdaPE = 8.65e-10
dM = 0.015
dP = 1.73
M0 = 5e-5
kaic = 1
KLox = 0.5



P1 = FiniteElement('P', tetrahedron, 1)
element = MixedElement([P1, P1, P1, P1, P1, P1])
Mixed_Space = FunctionSpace(Domain, element)
v1 , v2, v3, v4, v5, v6=TestFunction(Mixed_Space)
U = Function(Mixed_Space)
L, H, r, Lox, M, PP = split(U)

F=kl*r*L*v1*dx-Dl*alphaL*(L0-L)*v1*ds(1)+Dl*dot(grad(L),grad(v1))*dx\
 +kh*r*H*v2*dx-Dh*alphaH*(H0-H)*v2*ds(1)+Dh*dot(grad(H),grad(v2))*dx\
 +r*(kl*L+kh*H)*v3*dx-r0*v3*dx+Dr*dot(grad(r),grad(v3))*dx\
 +Lambda*M0*Lox*v4*dx-kl*r*L*v4*dx+Dlox*dot(grad(Lox),grad(v4))*dx\
 + dM*M*v5*dx-Dm*alphaM*M0*v5*ds(1)+Dm*alphaM*M*v5*ds(1)+Dm*dot(grad(M),grad(v5))*dx\
 +kaic*M*dot(grad(PP),grad(v5))*dx\
 + dP*PP*v6*dx-LambdaPE*(Lox/(KLox+Lox))*v6*dx+Dp*dot(grad(PP),grad(v6))*dx

solve(F == 0, U, solver_parameters={"newton_solver": {"relative_tolerance": 1e-12,   
       "absolute_tolerance": 1e-12,"maximum_iterations": 500}}) ​
and the mesh files are:
File attached: general_boundary_wenrui+5.xml (1.36 MB)

File attached: general_boundary_wenrui+5_facet_region.xml (1.1 MB)

File attached: general_boundary_wenrui+5_physical_region.xml (509.85 KB)

Community: FEniCS Project
You need to provide more code so we understand what your terms are. Also, your weak form seems wrong. What's the expression for M? You might find this example useful https://fenicsproject.org/pub/tutorial/sphinx1/._ftut1005.html#defining-subdomains-for-different-materials
written 9 weeks ago by Miguel  
Thanks for the reply. I added my code and I fixed that weak form (which is on the 6th line of my weak form), but I still receive that error.
written 9 weeks ago by Nima  
Regarding the weak form you have earlier mentioned, why don't you try  $M\nabla Pv.n|_{surface}-\int M\nabla P.\nabla vdx$MPv.n|surƒ aceMP.∇vdx.

Also try using nabla_grad(v) instead of grad(v) for the volume integral term.
written 9 weeks ago by Raghu Nandan  
1
About the surface term, it goes to zero because of my Neumann BC, But could you explain why I should use nabla_grad()? I was under the impression that it was for vector functions, but here all of my functions are scalar.
written 9 weeks ago by Nima  
Yeah. I am sorry. I haven't noticed that you are using scalars. grad() should work for you.

If you are still not able to solve the issue, try inner() instead of dot(). Essentially both do the same thing, but I use inner(). Maybe it will work for you too.
written 9 weeks ago by Raghu Nandan  
Is M a scalar or a tensor?
Why would this not be your form (Poisson's eqn., M scalar or tensor):
a = inner(grad(v), M*grad(p))*dx
written 9 weeks ago by jwinkle  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »