# Simulating advection - diffusion problem in a network of 1D pipe

I'm interested in solving the following advection-diffusion system in a 1D network of pipes.

$$\frac{\partial C}{\partial t} = D\frac{\partial ^2 C}{\partial x^2} - v\frac{\partial C}{\partial x}$$

To solve for the concentration of a solute in a 1D pipe of uniform dimension, without branches, I could use MATLAB's pdepe solver.

Could someone suggest if there is any code/tool that is available for solving advection/diffusion problem in a 1D network?

• I guess as soon as you have branches, you are in a quasi-2D setting. Is it possible to write it as a system of pde's where you have the branching points as the couplings between your equations? I'm not knowledgable in matlab, but I would presume that there is some functionality to solve a system of coupled equations. The tricky part will be how to express the couplings. In any branching point, the dynamics are not only influenced by the two neighboring points, but by a third. Do you have any sketch, or information on your pipe-network? – MPIchael Nov 5 '19 at 10:35
• you might also need some condition on your advection speed v. If the fluxes in one branching point don't add up, you will run into problems – MPIchael Nov 5 '19 at 10:37
• @MPIchael Please find the sketch of my network here. Yes, I could use the pdepe solver in MATLAB to solve the coupled equations. I am not really sure how to define flux conservation. I'd also like to know if there are solvers that can be used for this kind of problem. – Natasha Nov 5 '19 at 11:39
• In electronic circouit networks, you have certain rules about the fluxes (Kirchoff's). The fluxes in your case should behave similar, I think. One of the rules is, that at a node, all the incoming/outgoing fluxes have to add to zero. Otherwise you would have some sort of sink, or source in your network. – MPIchael Nov 5 '19 at 12:28

Based on the image that you provided in your comment, I believe you formulate your problem as a system of PDEs for each branch in your network and make sure at each connecting node mass is conserved. Let's say you have $$N$$ branches, so you need to solve the system of advection-diffusion equations for each branch ($$1 \leqslant i \leqslant N$$):

$$\frac{\partial C_{i}}{\partial t} + v_{i} \frac{\partial C_{i}}{\partial x} = D \frac{\partial C_{i}}{\partial x^{2}}$$

Let's say you have $$I$$ inflow points and $$O$$ outflow points. I can map these inflow and outflow points to IDs of branch by introducing $$F_{inflow}$$ and $$F_{outflow}$$ functions. In fact you set the ID of your inlet or outlet and $$F$$ will give you the ID of branch:

$$F_{inflow}(I_{v}) = i_{inlet}$$

$$F_{outflow}(O_{v}) =i_{outlet}$$

Where $$1 \leqslant I_{v} \leqslant I$$, $$1 \leqslant O_{v} \leqslant O$$, and $$i_{inlet}$$ and $$i_{outlet}$$ are IDs of branches that have inlet or outlet respectively.

I assume you put constant concentration at the inlets and zero diffusive flux at the outlets, but you can change it if you want:

$$C_{i_{inlet}}|_{x_{i_{inlet}}} = C_{0}^{i_{inlet}}$$

$$\frac{\partial C_{i_{outlet}}}{\partial x}|_{x_{i_{outlet}}} = 0$$

Where $$x_{i_{inlet}}$$ and $$x_{i_{outlet}}$$ are coordinates of $$i_{inlet}$$th and $$i_{outlet}$$th inlet and outlet respectively.

Now, let's say you have $$M$$ connecting points. At $$m$$th connecting point ($$1 \leqslant m \leqslant M$$), you have $$G(m)$$ branches. Conservation of mass should hold at each connecting point:

$$\sum_{\nu = 1}^{G(m)} J_{\nu} = 0$$

Where $$J_{\nu} = -D \frac{\partial C_{\nu}}{\partial x} + v_{\nu} C_{\nu} |_{x = x_{m}}$$ is the flux and $$x_{m}$$ is the coordinate of $$m$$th connecting point. Note that the flux despite its 1D nature has direction and as a result you don't need to distinguish incoming fluxes from outcoming fluxes in the above summation. Finally, you need to discretize this system of PDEs with their corresponding boundary conditions by something like finite difference in the simplest case. I don't think you can something ready in the MATLAB to feed this system of PDEs into it and then extract the results. So, it needs some efforts to implement this scheme and verify it with some simple analytical solution. I suggest to look at some electrical circuit solvers to get an idea but keep in mind due to the fact that you have diffusion here, you can’t map your problem into a pure electrical circuit and use existing electrical circuit solvers, unless you want to ignore diffusion.

• Just to clarify, do electrical circuit solvers consider diffusion? – MPIchael Nov 6 '19 at 13:21
• @MPIchael No, sorry that was a typo... – Alone Programmer Nov 6 '19 at 13:25
• @AloneProgrammer When we say $J_v$ is the flux, could you please clarify whether $J_v$ is the flux that is entering or leaving a node? – Natasha Nov 16 '19 at 12:08
• @Natasha "Note that the flux despite its 1D nature has direction and as a result you don't need to distinguish incoming fluxes from outcoming fluxes in the above summation." – Alone Programmer Nov 16 '19 at 18:16
• @AloneProgrammer Thank you for the response. When we say the summation of flux = 0 at the connecting points, do we assume connecting points have negligible volume? Also, should$C_v$ be written in the units of quantity/length? quantity here refers to mass/mole. Sorry, I am confused about the units of flux. Is it mole/time? – Natasha Nov 17 '19 at 2:09