Reaction-advection-diffusion in VoronoiFVM

General Usage
1

I was trying to use VoronoiFVM.jl and I am not sure I understand how to build a Physics object correctly.
This is probably because of my lack of knowledge about FVM in general.

I have a PDE in 2D with two species and the following shape:

\begin{aligned} \nabla_t u &= \nabla_x (D_u \nabla_x u - k_1 u \nabla_x v) + (k_2 (1 - u) - k_3 + f(v)) u, \\ \nabla_t v &= D_v \nabla_x v + k_4 u - k_5 v + k_6. \end{aligned}

Would it be correct to translate it to this Physics object?

grid = simplexgrid(collect(0:1/100:1), collect(0:1/100:1))

physics = VoronoiFVM.Physics(
                         ; reaction = function (f, u, node, data)
                             f[1] = (k_2*(1 - u[1]) - k_3 + f(u[2]))*u[1]
                             f[2] = k_4*u[1] - k_5*u[2]
                         end,
                         flux = function (f, u, edge, data)
                             nspecies = 2
                             f[1] = D_u*(u[1,1]-u[1,2]) - k_1*u[1]*(u[2,1]-u[2,2])
                             f[2] = D_v*(u[2,1]-u[2,2])
                         end, source = function (f, node, data)
                             f[1] = 0
                             f[2] = k_6
                         end, storage = function (f, u, node, data)
                             f[1] = u[1]
                             f[2] = u[2]
                         end)

sys = VoronoiFVM.System(grid, physics; unknown_storage = unknown_storage)

enable_species!(sys, 1, [1])
enable_species!(sys, 2, [1])

boundary_neumann!(sys, 1, 1, 1.0)
boundary_neumann!(sys, 1, 2, 0.0)
boundary_neumann!(sys, 2, 1, 1.0)
boundary_neumann!(sys, 2, 2, 0.0)
1 Like
3

Hi,
this already looks quite good.

Two things:

  • VoronoiFVM assumes that all u-dependent terms are on the left side of the equal sign. So you should reverse the sign in the reaction term.
  • For the convective term, you should try some upwinding. The first equation looks similar to a drift-diffusion equation where the drift (convection) term involved a gradient. The best way to discretize this is to use an exponential fitting ansatz. You might consider to look at classflux in 160: Unipolar degenerate drift-diffusion · VoronoiFVM.jl
4

Sorry for the late answer.

Yes, the reactions have the sign reversed. Thanks for pointing out that.

I’ve tried reading more about FVM, but I do not fully understand how it would be implemented in this case.

From what I see in the example, classflux! seems to only be defined and does not appear anywhere else.
But what I see from that and 103: 1D Convection-diffusion equation · VoronoiFVM.jl, is that a bernoulli function is used. Is that necessary?

I have this, but the advection term seems to not be correct:

physics = VoronoiFVM.Physics(
                         ; reaction = function (f, u, node, data)
                             f[1] = k_3 - (k_2*(1 - u[1]) - f(u[2]))*u[1]
                             f[2] = k_5*u[2] - k_4*u[1]
                         end,
                         flux = function (f, u, edge, data)
                             nspecies = 2
                             f[1] = D_u*(u[1,1]-u[1,2]) - k_1*u[1]*(u[2,1]-u[2,2])
                             f[2] = D_v*(u[2,1]-u[2,2])
                         end, source = function (f, node, data)
                             f[1] = 0
                             f[2] = k_6
                         end, storage = function (f, u, node, data)
                             f[1] = u[1]
                             f[2] = u[2]
                         end)

Is there a way to solve it? I am interested in temporal dynamics.