# How to solve advective equation with source term depending on variable

I have the following equation $$\dfrac{\partial s}{\partial t} + \nabla \cdot \left( \vec{v} s\right) = f(s)$$

Where $$f(s)$$ is an explicit source term that depends on $$s$$, e.g., ($$\sin(s)\;cos(s)$$).

When solving with finite volume method (fvm), should I solve the above equation once until a final tolerance? Or should I keep on solving the equation until the initial residual is below a tolence?

Basically, should it be:

1:

Begin time step

Assemble linear system of Eqs and solve for $$s$$

End time step


or 2:

Begin time step

while(EqInitialResidual > tolerance)
Assemble linear system of Eqs and solve for $$s$$
(Here, solve once for $$s$$. Take the new value of $$s$$, compute $$f(s)$$, solve again for $$s$$ until the result is bellow a tolerance.
End time step


What is the best solution approach?

• Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking.
– Community Bot
Commented Sep 8, 2022 at 16:49

It depends on your temporal discretization. If you have an explicit time discretization, then $$s_{n+1} = \mathcal{L}(s_n)$$ Since $$s_n$$ is known, you can explicitly compute $$s_{n+1}$$ using the operator $$\mathcal{L}(s)$$, and this can be evaluated in a bounded finite number of steps.

For example, suppose you only have a 1D problem and are using the forward Euler method and an upwinding flux (assuming that $$v > 0$$). Then $$s_{n+1}(x_i) = s_n(x_i) + \Delta t \left(f(s_n(x_i)) + \frac{v}{\Delta x} (s_{n}(x_{i-1}) - s_{n}(x_{i}))\right)$$ Simply plug in $$s_n$$ to quickly compute $$s_{n+1}$$. Trying to iterate this multiple times will not improve your solution since you will get the same result every time.

If you have an implicit temporal discretization, this is no longer true as we now have $$\mathcal{L}(s_n, s_{n+1}) = 0$$ where generally there is no way to express $$s_{n+1}$$ purely in terms of $$s_n$$. The most common solution to this is to use the Newton-Raphson method to break your problem up into a bunch of quasi-linear problems, and iterate until you find a solution $$u_{n+1}$$ which sufficiently minimizes the residual (what you described as method 2).

Repeating our previous 1D example, but this time using Backwards Euler: $$\mathcal{L}(s_n,s_{n+1}) = 0 = s_{n+1}(x_i) - s_n(x_i) - \Delta t \left(f(s_{n+1}(x_i)) + \frac{v}{\Delta x} (s_{n+1}(x_{i-1}) - s_{n+1}(x_{i}))\right)$$ This time it's not possible to explicitly solve for $$s_{n+1}$$ in terms of $$s_n$$ and constants (or if it is possible it will be relatively difficult), so we must rely on some non-linear root solving algorithm, which often iterate until a solution sufficiently close to the true solution is found.

For the Newton-Raphson method, this looks like: $$\partial_{s'_{n+1}(x_i)}\left(\mathcal{L}(s_n,s'_{n+1})\right) \Delta s'_{n+1}(x_i) = s'_{n+1}(x_i) - s_n(x_i) - \Delta t \left(f(s'_{n+1}(x_i)) + \frac{v}{\Delta x} (s'_{n+1}(x_{i-1}) - s'_{n+1}(x_{i}))\right)$$ where $$s'$$ is the current guess for $$s_{n+1}$$, and the next guess for $$s_{n+1}$$ is $$s''_{n+1} = s'_{n+1} + \Delta s'_{n+1}$$ Repeat this until $$\|\mathcal{L}(s_n, s'_{n+1})\| < \epsilon$$. This could take one iteration, or it may take hundreds.

Side-note: This process could be non-convergent, though explaining methods for improving the Newton-Raphson method to achieve global convergence is beyond the scope of this answer. Search for "line search" or "trust regions" in the context of nonlinear optimization for more information on this.

• I assume the trade off between the approaches is the calculation stability. Although the explicit method will compute things in one step it will be unstable (since it should require a smaller Delta t.
– mysn
Commented Sep 8, 2022 at 8:28
• Yes, all explicit methods have a maximum stable timestep, while some implicit methods theoretically are stable for any $\Delta t$. However, this doesn't necessarily cover everything since as you take larger $\Delta t$ solving the implicit system becomes harder, especially for non-linear problems. Both methods are also subject to truncation errors which scale with $\Delta t$, so you might need to take a smaller $\Delta t$ to get the desired accuracy. Commented Sep 8, 2022 at 8:32
• In this case I am using implicit formulation for the time-rate-of-change and for the divergence and explicit calculation of the f(s). In the end it will result in a "mixed" type formulation. However, I have been using the 2nd method to solve the system to 1e-5 which I do not know if it is good enough, or if I should switch to full explicit formulation and solve things in one step.
– mysn
Commented Sep 8, 2022 at 8:55
• Depending on your situation I would highly recommend looking at operator splitting approaches. They would allow you to split the problem up into explicit time stepping for the advection dominated problem and then utilize an implicit method for the source term in a very straight forward manner. There are also IMEX methods that do this in a more sophisticated way as well. Commented Sep 8, 2022 at 15:34
• Whether fully implicit or mixed implicit/explicit would be beneficial depends on the characteristics of $f$ vs. your desired spatial resolution/value of $\vec{v}$. If $\partial_s f \approx \frac{v \Delta t}{\Delta x}$, then you likely won't benefit much from using implicit methods. Commented Sep 8, 2022 at 20:43