Solving a PDE with a global (integral) constraint

[Amit Singh]

Hil

I am trying to solve a reaction-diffusion PDE with two reactants A and B. Reactant B diffuses much faster than A i.e. it instantly diffuses to become uniform over the whole domain in response to changes in concentration of A. Thus, the concentration of this reactant is simply a real number that varies with time but not with space. The total mass of reactants A and B is fixed. The governing equations are as shown in the following image

Here, a(x, t) and b(t) are the concentrations of the reactants. D is the diffusion co-efficient and f(a, b) is a non-linear reaction term. M is total mass of A and B and S is the surface area of the domain. I came up with the following weak form from these two equations

In this weak form, `u` and `r` are the test functions corresponding to `a` and `b`. I have solved this problem in FEniCS FEM software already but I need to implement adaptive mesh refinement for which I prefer deal.ii. In FEniCS, I created a mixed finite element space consisting of P1 and R finite element spaces. R is the simply the space of real numbers which agrees with `b` being a single value for the whole domain at time `t`. I have been reading the deal.ii tutorials and documentations. I have gathered that I will need to use FESystem to compose two finite element spaces for `a` and `b`. I can use FE_Q(1) for `a` but I am not able to find what I should use for `b`.  Please help.

Thanks and regards,
Amit

[Wolfgang Bangerth]

Amit,

> I am trying to solve a reaction-diffusion PDE with two reactants A and B.
> Reactant B diffuses much faster than A i.e. it instantly diffuses to become
> uniform over the whole domain in response to changes in concentration of A.
> Thus, the concentration of this reactant is simply a real number that varies
> with time but not with space. The total mass of reactants A and B is fixed.
> The governing equations are as shown in the following image
>

> System.png
>
> Here, a(x, t) and b(t) are the concentrations of the reactants. D is the
> diffusion co-efficient and f(a, b) is a non-linear reaction term. M is total
> mass of A and B and S is the surface area of the domain. I came up with the
> following weak form from these two equations
>

> WeakForm.png

You can do it that way, but (i) that's awkward because now you have this
global term in your problem, and (ii) it's not really useful to think of a
scalar equation as one you wanted to solve with a "field" for 'b' and 'r'.

But let's pretend that you wanted to do it that way (and that for simplicity,
f(a,b) is a linear function of the two variables, so that your problem is
linear). Then the way I would approach this is by conceptually saying that you
have one DoFHandler for 'a', which is your field, and your overall linear
problem leads to a matrix that you can partition as follows:

[ A C ] [ a ] = [ ... ]
[ B^T I ] [ b ] = [ M/S ]

Here, C results from f(a,b) = X a + C b, and B^T is a 1xn vector of all ones.
I is the 1x1 unit matrix. A is what you get from the spatial discretization of
'a' plus the X that comes out of 'f'.

You can build this linear system as a BlockSparseMatrix, where only A is
associated with the DoFHandler.


I would imagine that a better way to solve this problem is to do some sort of
operator splitting where you alternate between updating 'a' and 'b'. The
equations for the updates are obvious.

Best
W.

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/

[Amit Singh]

Thank you very much for the response Dr. Bangerth. I admire that you have been a lead developer for deal.ii for all this time and are still actively maintaining it.

I am new to operator splitting. I will try to read up about it.

[Wolfgang Bangerth]

On 5/6/20 9:47 AM, Amit Singh wrote:
> Thank you very much for the response Dr. Bangerth. I admire that you have been
> a lead developer for deal.ii for all this time and are still actively
> maintaining it.
> I am new to operator splitting. I will try to read up about it.

Thanks for the kind words.

Think about the approach I outlined like the one we use in step-21.