Как вы задаете граничное условие Неймана (фиксированный поток, нормальный к поверхности) в Fipy?

Как мне явно установить поток, нормальный к граничной грани в гибкой сетке, как определенное значение, не ограничивая компоненты потока внутри грани?

Граничное условие Неймана может быть определено как: (1) фиксированная составляющая потока, нормальная к граничной грани, или (2) как полная спецификация потока на грани. Условие fipy по умолчанию - первое (значение = 0), а явный метод (faceGrad.constrain) - второе. Эту проблему можно понять, добавив следующий код в конец файла fipy diffusion.mesh20x20 пример и отмечая разные результаты градиента лица.

facesNeumann = mesh.exteriorFaces & ~facesTopLeft & ~facesBottomRight
print 'grad(phi) with default Neumann BC...'
print phi.faceGrad.value.T[facesNeumann.value]
phi.faceGrad.constrain(0.,facesNeumann)
DiffusionTerm().solve(var=phi)
print 'and with explicit Neumann BC...'
print phi.faceGrad.value.T[facesNeumann.value]

python mesh pde fipy
person Rupert    schedule 18.09.2016    source источник

Ответы (3)


arrow_upward
2
arrow_downward

См. Обсуждение фиксированных граничных условий потока. По сути, вы добавляете источник, содержащий расхождение желаемого граничного потока, к условию отсутствия потока по умолчанию FiPy.

person jeguyer    schedule 19.09.2016

arrow_upward
2
arrow_downward

С точки зрения простого указания потока, нормального к границе, для переменной, не зависящей от решения уравнения, похоже, нет никакого способа сделать это. Например, синтаксис может быть phi.faceGrad[0].constrain(...), но в настоящее время он не работает в FiPy. Это также может быть сложно реализовать для произвольно ориентированных лиц.

Однако для практических целей компонент, касательный к границе, не используется при решении уравнения в FiPy, только нормальный компонент имеет какое-либо влияние на решение. Например, возьмите следующий код,

import fipy as fp
mesh = fp.Grid2D(nx=2, ny=1)
var = fp.CellVariable(mesh=mesh)
var.constrain(1, mesh.facesLeft)
var.constrain(0, mesh.facesRight)
#var.faceGrad.constrain(0, mesh.facesTop)
fp.DiffusionTerm().solve(var)
print 'face gradient on top plane:',var.faceGrad[0, mesh.facesTop.value]
print 'variable value:',var

Это дает результат

face gradient on top plane: [-0.5 -0.5]
variable value: [ 0.75  0.25]

Ответ правильный, но градиент верхней грани -0,5. Однако, когда строка #var.faceGrad.constrain(0, mesh.facesTop) не закомментирована, результат будет

face gradient on top plane: [ 0.  0.]
variable value: [ 0.75  0.25]

Тангенциальный градиент грани теперь равен 0, но ответ тот же. Дело в том, что установка градиента грани по касательной (через .constrain) не влияет на решение.

person wd15    schedule 19.09.2016
comment
Что делать, если вы хотите, чтобы одно граничное условие зависело от величины градиента другого, и поэтому вы хотите установить только нормальный компонент. Скажем, phi1.constrain(fp.numerix.dot(phi2.faceGrad, phi2.faceGrad)), где мы хотим указать только нормальный компонент phi2 на границе. Означает ли этот ответ, что это не сработает? (У меня нет хорошего примера, просто любопытно). Может быть, в таком случае уместен ответ @jguyer? - person muon; 20.09.2016
comment
@muon; Я думаю, что метод ограничения не работает для большинства произвольных сложных взаимозависимых граничных условий, и требуется метод исходных условий (ответ @jeguyer). - person wd15; 26.09.2016
comment
Имеет смысл. Спасибо! - person muon; 26.09.2016

arrow_upward
1
arrow_downward

Обсуждение граничных условий с фиксированным потоком действительно отвечает на мой вопрос, но я не понял документацию, которая очень краткая. Ниже приведен рабочий пример, который показывает, как применить это в простом двухмерном случае, похожем на пример с сеткой 20x20.

import matplotlib.pyplot as plt
from fipy import *

nx = 20
ny = nx
dx = 1.
dy = dx
L = dx * nx
msh = Grid2D(dx=dx, dy=dy, nx=nx, ny=ny)

xFace, yFace = msh.faceCenters
xCell,yCell = msh.cellCenters    
phi = CellVariable(name = "solution variable",
                   mesh = msh,
                   value = 0.)    
D = FaceVariable(name='diffusion coefficient',mesh=msh,value=1.)

# Dirichlet condition on top face
valueTop = FaceVariable(name='valueTop',mesh=msh,value=xFace*0.1-1)
phi.constrain(valueTop, msh.facesTop)

# Fixed flux (Neumann) on base
D.constrain(0., msh.facesBottom)
fluxBottom = -0.05

eq = DiffusionTerm(coeff=D) + (msh.facesBottom * fluxBottom).divergence
eq.solve(var=phi)

# confirm only y component of gradient is zero
print phi.faceGrad.value.T[msh.facesBottom.value]

plt.scatter(phi.value, yCell)
plt.show()

viewer = Viewer(vars=phi, datamin=-1., datamax=1.)
viewer.plot()
person Rupert    schedule 26.09.2016