jorgensd / jorgensd/MultiMeshShapeOpt_code

MultiMesh in solving eigenvalue problem of Stokes operator

オープン
#1 コメント 0 件 リアクション 0 件 担当者 0 名 GitHub で見る
主要言語
Python
スター
3
フォーク
1
PR マージ指標
30日以内にマージされた PR はありません

説明

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.)

コントリビューションガイド

このリポジトリのコントリビューションガイドは索引されていません

評価

この issue はまだ評価されていません。

新しい issue をメールで受け取る

初心者向けの GitHub issue を短くまとめたダイジェスト。