On Aug 2, 9:12Â pm, Gordon Sande
worldnet.att.net> wrote:> On 2008-08-02 15:32:50 -0300, s.c.kra...@
gmail.com said:
>
>
>
>> On Aug 2, 5:21Â pm, boulou...@
gmail.com wrote:
>>> On 31 juil, 20:49, boulou...@
gmail.com wrote:
>
>>>> I am working on a 3d finite volume scheme for an advection-diffusion-
>>>> reaction problem involving a large number of chemical species (more
>>>> than 60) and a large domain (an big lake for example). Since this
>>>> scheme will be used on large problem, I want it to be as efficient as
>>>> possible. The linear operators are splitted in 2 :
>
>>>> (1) advection-diffusion is solved using a fully implicit finite volume
>>>> discretisation with a multigrid method for solving the linear system
>>>> of equations
>>>> (2) chemistry is solved using a Runge-Kutta-Rosenbrock solver for
>>>> stiff ODE.
>
>>>> The transport (1) actually have the following form
>
>>>> do specie=1, nbSpecies
>>>> Â Â Â Â call construct_matrix(specie)
>>>> Â Â Â Â call solve_linear_system(specie)
>>>> end do
>
>>>> and takes a lot of time on the computer.
>
>>>> Assuming that diffusion coefficients are the same for all species, the
>>>> whole fluid (including all species) should follow the same path during
>>>> the transport. I wonder if it really necessairy to loop over all
>>>> species and compute the transport several time. It is possible to
>>>> compute the transport of the fluid once, and after reuse this
>>>> calculation to the different species ?
>
>>>> I would really appreciate suggestion or reference on this.
>
>>> Maybe I could do the following :
>
>>> 1) Suppose I use a fully implicit (backward Euler for example) or semi-
>>> implicit (Crank-nicholson) time integrating scheme. I will use a
>>> finite volume discretisation to create a linear system in the form
>
>>> A * c(n+1) = c(n)
>
>>> 2) Compute
>
>>> T = inverse(A)
>
>>> Using a finite volume discretisation, with an hybrid differencing
>>> scheme (central/upwind depending on the Peclet number) for the
>>> advective part, the matrix A will be diagonnally dominant and all
>>> entries are positive. Hence, the inverse of A exists and could be
>>> computer with an appropriate method (any suggestions on an efficient
>>> way to do this ?)
>
>>> 3) loop over all species and solve
>
>>> C(n+1) = T * C(n)
>
>>> for each one.
>
>> I think you're missing an important point here. First of all inverting
>> a matrix gives you a dense matrix, as opposed to your matrix A which
>> will be sparse. Even for a very small number of nodes/cells/matrix
>> rows of 10.000 you already need 800MB to store it. More importantly
>> inverting matrices is computationally expensive, where a single matrix
>> equation can be solved in O(n) to O(n^2) time, inverting the whole
>> matrix costs O(n^2) to O(n^3). So this is only worth it if the number
>> of rhs's is in the order of the number of rows of your matrix.
>
>> Cheers
>> Stephan
>
> Two points:
>
> 1. If A is sparse then A^-1 will almost surely be dense BUT A=LU may lead
> Â Â to lower density factors. Solving with LU is as easy as multipying by
> Â Â A^-1 and finding LU is cheaper that finding A^-1. Nice sparsity structure
> Â Â and a good sparse solver can often do surprising well as solutions tend
> Â Â to be local. This is a fact that can be hard to exploit other ways. Details
> Â Â depend on particular problems.
>
Heh yes, I was oversimplifying things and merely responding to the
OP's suggestion of computing the explicit inverse of his matrix.
Pointing at LU would of course have been much better advice.
Stephan
> 2. You should not be looking for numerical analysis advice in a programming
> Â Â newsgroup. c.l.f is the least bad that way as many of the regulars are
> Â Â knowledgable in numerical analysis.
>
> Both points are extremely standard advice.
