Galerkin / Poisson / Fenics Tidak Terputus

10

Saya mencoba menyelesaikan persamaan 2D Poisson menggunakan metode Discontinuous Galerkin (DG) dan diskritisasi berikut (saya punya file png tapi saya tidak diizinkan mengunggahnya, maaf):

Persamaan:

(κT)+f=0

Persamaan baru:

q=κTq=f

Bentuk lemah dengan numerik fluks T dan q :T^q^

qwdV=T(κw)dV+κT^nwdSqvdV=vfdV+q^nvdS

Fluks numerik (metode IP):

q^={T}C11[T]T^={T}

dengan

{T}=0.5(T++T)[T]=T+n++Tn

Saya menulis skrip python fenics sederhana untuk menyelesaikan persamaan. Solusi yang saya dapatkan tidak baik. Saya akan sangat menghargai jika seseorang yang akrab dengan metode DG bisa melihat skrip di bawah ini dan memberi tahu saya apa yang saya lakukan salah.

Formulasi galerkin terus menerus yang saya tambahkan dalam skrip memberikan solusi yang bagus.

Terima kasih banyak sebelumnya.

from dolfin import *

method = "DG" # CG / DG

# Create mesh and define function space
mesh = UnitSquare(32, 32)
V_q = VectorFunctionSpace(mesh, method, 2)
V_T = FunctionSpace (mesh, method, 1)
W = V_q * V_T

# Define test and trial functions
(q, T) = TrialFunctions(W)
(w, v) = TestFunctions(W)

# Define mehs quantities: normal component, mesh size
n = FacetNormal(mesh)

# define right-hand side
f = Expression("500.0*exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)")

# Define parameters
kappa = 1.0

# Define variational problem
if method == 'CG':
  a = dot(q,w)*dx \
       + T*div(kappa*w)*dx \
       + div(q)*v*dx

elif method == 'DG':
  #modele = "IP"
  C11 = 1.

  a = dot(q,w)*dx + T*div(kappa*w)*dx \
      - kappa*avg(T)*dot(n('-'),w('-'))*dS \
                                           \
      + dot(q,grad(v))*dx \
      - dot( avg(grad(T)) - C11 * jump(T,n) ,n('-'))*v('-')*dS

L = -v*f*dx

# Compute solution
qT = Function(W)
solve(a == L, qT)

# Project solution to piecewise linears
(q , T) = qT.split()

# Save solution to file
file = File("poisson.pvd")
file << T

# Plot solution
plot(T); plot(q)
interactive()
micdup
sumber

Jawaban:

10

Untuk menerapkan masalah Anda dalam FEniCS, Anda harus mengganti integral dalam hal batas dengan integral dalam hal tepi. Ini memperkenalkan lompatan / rata-rata dalam fungsi tes, yang Anda sepenuhnya lewatkan dalam implementasi Anda. Oleh karena itu, sistem tidak dapat dibalik dan solusi Anda tidak beres. Persamaan (3.3) dalam Arnold et. Al. 2002 memberi Anda alat untuk menulis ulang formulir yang lemah:

KThKqKnKϕKds=Γ[q]{ϕ}ds+Γ0{q}[ϕ]ds

Di sini adalah penyatuan tepi Anda dan sama tanpa batas.ΓΓ0

Sekarang fluks Anda bernilai tunggal, yang berarti Anda dapat menjatuhkan lompatan fluks Anda. Karenanya

KThKq^nKvKds=Γ0q^[v]ds+Ωq^nvdsKThKwnKκT^ds=Γ[w]κT^ds

Ini membawa kami ke modifikasi kode Anda berikut:

C11 = 1.
qhat = avg(grad(T)) - C11 * kappa*jump(T,n)
qhatbnd = grad(T) - C11 * kappa*T*n

a = dot(q,w)*dx + T*div(kappa*w)*dx \
  - kappa*avg(T)*jump(w,n)*dS \
  - kappa*T*dot(w,n)*ds \
  - dot(q,grad(v))*dx \
  + dot( qhat, jump(v,n))*dS \
  + dot( qhatbnd, v*n)*ds

Saya belum punya waktu untuk benar-benar mencoba ini, jadi berhati-hatilah dengan kemungkinan kesalahan tanda-tanda dll. Tapi saya harap ini bisa membantu.

Referensi: DN Arnold, F. Brezzi, B. Cockburn, LD Marini: Analisis Terpadu Metode Galerkin Terputus untuk Masalah Elliptic SIAM J. Num. Anal, 39 (2002), 1749-1779

Christian Waluga
sumber
Ya saya benar-benar kehilangan sesuatu.
micdup
-2

Ya saya benar-benar kehilangan sesuatu!

Sekarang berfungsi dengan baik.

Terima kasih banyak atas bantuannya!

micdup
sumber
2
Demi kelengkapan, bisakah Anda menggambarkan apa yang hilang dan bagaimana Anda memperbaikinya.
Paul