### solve Anisotropy equation

197

views

-2

Hello everyone!

∇.(σ∇u)=0

In this problem, σ Is a tensor （σ is as_matrix([[1, 2],[2, 3]])）

How should I write this variational equation and how to solve it？

thank you!

∇.(σ∇u)=0

In this problem, σ Is a tensor （σ is as_matrix([[1, 2],[2, 3]])）

How should I write this variational equation and how to solve it？

thank you!

Community: FEniCS Project

### 1 Answer

7

You can modify the Poisson demo here

http://fenics.readthedocs.io/projects/dolfin/en/stable/demos/poisson/python/demo_poisson.py.html

by changing the bilinear form to

However, note that \(\sigma\) should be positive definite for the problem to be stable. (Physically: "heat needs to flow from hot to cold".) The suggested matrix has a negative eigenvalue, and will produce a very strange-looking discrete solution.

http://fenics.readthedocs.io/projects/dolfin/en/stable/demos/poisson/python/demo_poisson.py.html

by changing the bilinear form to

`a = inner(sigma*grad(u), grad(v))*dx`

However, note that \(\sigma\) should be positive definite for the problem to be stable. (Physically: "heat needs to flow from hot to cold".) The suggested matrix has a negative eigenvalue, and will produce a very strange-looking discrete solution.

1

You did not incorporate what David Kamensky wrote. The expression inside your bilinear form:

`inner(tao,grad(u))`

does not make sense. There is no inner product between tensors and vectors. Use `inner(tao*grad(u), grad(v))*dx`

written
4 months ago by
klunkean

By the modification, the following form:

inner(tao*grad(u), grad(v))*dx

L = f*v*dx + g*v*ds

but the running is error

inner(tao*grad(u), grad(v))*dx

L = f*v*dx + g*v*ds

but the running is error

written
4 months ago by
kevin

What error is the running? Maybe paste the message?

written
4 months ago by
klunkean

Use

`tao = as_tensor([[1,2],[3,4]])`

written
4 months ago by
Adam Janecka

thanks a lot! successful

written
4 months ago by
kevin

Please login to add an answer/comment or follow this question.

the demo is following:

from dolfin import *

import numpy as np

import matplotlib.pyplot as plt

mesh = UnitSquareMesh(32, 32)

V = FunctionSpace(mesh, "Lagrange", 1)

def boundary(x):

return x[0] < DOLFIN_EPS or x[0] > 1.0 - DOLFIN_EPS

u0 = Constant(0.0)

bc = DirichletBC(V, u0, boundary)

u = TrialFunction(V)

v = TestFunction(V)

f = Expression("10*exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)", degree=2)

g = Expression("sin(5*x[0])", degree=2)

tao = np.array([[1,2],[3,4]])

a = inner(inner(tao,grad(u)), grad(v))*dx

L = f*v*dx + g*v*ds

u = Function(V)

solve(a == L, u, bc)

plot(u)

plt.show()

interactive()

how can i solve the equation?