# Convection-diffusion problem with pre-computed velocity field

**URL:** <https://community.freefem.org/t/convection-diffusion-problem-with-pre-computed-velocity-field/630>\
**Category:** General Discussion\
**Created:** [November 1, 2020, 8:32pm UTC](https://community.freefem.org/t/convection-diffusion-problem-with-pre-computed-velocity-field/630 "2020-11-01T20:32:55Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![G\_B](https://avatars.discourse-cdn.com/v4/letter/g/779978/32.png) [@G\_B](https://community.freefem.org/u/G_B)\
**Post date:** [November 1, 2020, 8:32pm UTC](https://community.freefem.org/t/convection-diffusion-problem-with-pre-computed-velocity-field/630/1 "2020-11-01T20:32:55Z")

</div>

Dear FreeFem users,

I’m trying to solve a convection-diffusion problem of this form:

\partial\_t c - K\triangle c + \nabla\cdot(u c) - B(cv - c)=0

where K is the diffusion coefficient and it is constant, the term B(cv - c) is a distributed source in a region of the domain and u is the velocity field which has been previously computed (from Darcy’s law). In figure you can see u.

 ![Screenshot 2020-11-01 at 21.28.55](https://canada1.discourse-cdn.com/flex030/uploads/freefem/original/1X/9054bc058341d01bbb5479c8f4845e72d42d62ff.png)

I have tried to solve the equation in this form:  
\partial\_t c - K\triangle c + u\cdot\nabla c + (\nabla\cdot u )c - B(cv - c)=0

by using the Freefem convect operator on the time derivative and the term u\cdot\nabla c.  
Here it is the piece of code:

macro grad(fun) [dx(fun), dy(fun)] //eom  
Vh divvel = dx(uxh) + dy(uyh);  
problem Cevol(C, v, solver=GMRES) =  
int2d(Th)(1./dt_C_v + K_grad©’grad(v))  
- int2d(Th)(1./dtconvect([uxh, uyh], -dt, Cold)v)  
- int2d(Th)(Colddivvel_v)  
- int2d(Th)(B_Cv_v)  
+ int2d(Th)(B_C_v)  
;  
I would like to know if this is the correct way to treat the convection term of such an equation. The thing is that in my results the problem is solved with an acceptable solution when I discard the term - int2d(Th)(Cold_divvel_v), which is not correct in my opinion, since I do not have the condition \nabla\cdot u =0.  
If I consider the - int2d(Th)(Cold_divvel_v) in the problem, the solution explodes in few iterations in time.

Do you have any suggestion?

Thank you very much.

Best regards,  
Gabriele

---

<div class="post-metadata">

**Author:** ![frederichecht](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/frederichecht/32/15_2.png) [@frederichecht](https://community.freefem.org/u/frederichecht)\
**Post date:** [November 3, 2020, 10:52am UTC](https://community.freefem.org/t/convection-diffusion-problem-with-pre-computed-velocity-field/630/2 "2020-11-03T10:52:23Z")

</div>

In this schema explicite in the part c(\nabla.u) if (\nabla.u) is large the the schema  
become instable son the first thing is to implicite this term and warning to with the signe of  
int2d(Th)(Cold_divvel_ v)  
it is a + not a -.

```
problem Cevol(C, v, solver=sparsesolver) =
int2d(Th)(1./dt*C* v + K*grad(c)'*grad(v) + C*divvel* v)
- int2d(Th)(1./dt*convect([uxh, uyh], -dt, Cold)*v)
- int2d(Th)(B*Cv* v)
+ int2d(Th)(B*C* v)
;

```

Best Regards, F. Hecht