jorgensd / jorgensd/MultiMeshShapeOpt_code

MultiMesh in solving eigenvalue problem of Stokes operator

Abierto
#1 0 comentarios 0 reacciones 0 asignados Ver en GitHub
Lenguaje dominante
Python
Estrellas
3
Forks
1
Métricas de merge de PR
Sin PR fusionados en 30 d

Descripción

I am following your code to use multimesh to solve eigenvalue problem of the Stokes operator.

By selecting test function from a rawer mesh, I want to weaken the divergence-free condition.
When ``mesh_u`` and ``mesh_div_free`` are the same mesh, the eigenfunction should be completely divergence-free.
But the code below failed to apply the divergence-free condition to solution ``u``.

Could you give me any hint where the code is wrong?

```
from fenics import *

mesh_u=UnitSquareMesh(32, 32, diagonal="right")
mesh_div_free=UnitSquareMesh(16, 16, diagonal="right")
mesh_R=UnitSquareMesh(1, 1, diagonal="right")

multimesh = MultiMesh()
multimesh.add(mesh_u)
multimesh.add(mesh_div_free)
multimesh.add(mesh_R)
multimesh.build()

degree=3

CG = VectorElement("CG", triangle, degree,2)
DG = FiniteElement("DG", triangle, degree-1)
Real = FiniteElement("Real", triangle, 0)
MixedV = MultiMeshFunctionSpace(multimesh, MixedElement([CG, DG, Real]))
V = MultiMeshSubSpace(MixedV, 0)

(u, p, alpha) = TrialFunctions(MixedV)
(v, q, beta) = TestFunctions(MixedV)

a = inner(grad(u), grad(v))*dx + div(v)*p*dx + q*div(u)*dx + alpha*q*dx + beta*p*dx
b = inner(u,v)*dx

f = Constant((0.0, 0.0))
L = inner(f,v)*dx

mvc = MeshValueCollection("size_t", mesh_u, 1)
mf = cpp.mesh.MeshFunctionSizet(mesh_u, mvc)

noslip = Constant((0.0, 0.0))
marker=100
bc0=MultiMeshDirichletBC(V, noslip, mf, marker, 0)

A = assemble_multimesh(a)
B = assemble_multimesh(b)

bc0.apply(A)
bc0.zero(B)

# downcast to PETSc matrices
MatA = as_backend_type(A)
MatB = as_backend_type(B)

import scipy.io as myio
import scipy.sparse as sp

row_, col_, val_ = MatA.mat().getValuesCSR()
mat_A = sp.csr_matrix((val_,col_,row_))
row_,col_,val_ = MatB.mat().getValuesCSR()
mat_B = sp.csr_matrix((val_,col_,row_))

#Another approach to remove boundary value DOF
bd_int = inner(u,v)*ds
bd_int = assemble_multimesh(bd_int)
Mat_bd_int = as_backend_type(bd_int)

row_,col_,val_ = Mat_bd_int.mat().getValuesCSR()
Mat_bd_int = sp.csr_matrix((val_,col_,row_))
myio.savemat("./matrix/matrix.mat",{'A':mat_A,'B':mat_B,'BD_Int':Mat_bd_int})
```
(I exported the matrix to mat file and then solve it in Octave or MATLAB.)

Guía de contribución

No hay ninguna guía de contribución indexada para este repositorio

Evaluación

Este issue todavía no se ha evaluado.

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.