-------- Forwarded Message --------
Subject: Re: [OPM] Making OPM modular: support for more linear solvers
Date: Tue, 28 Oct 2014 12:35:26 +0000
From: Atgeirr Rasmussen <[email protected]>
To: Jørgen Kvalsvik <[email protected]>
28. okt. 2014 kl. 13:27 skrev Jørgen Kvalsvik <[email protected]>:
On 10/28/2014 01:19 PM, Atgeirr Rasmussen wrote:
28. okt. 2014 kl. 12:14 skrev Jørgen Kvalsvik <[email protected]>:
Hi,
I am a student at NTNU http://www.ntnu.no/, and my pre-master project is
modularizing OPM by implementing (at least trying to) support for using more
linear solvers. There is some code available, see the LinearSolverInterface and
friends in Core, but there is very little code that actually uses this. The
biggest sinner I've found so far is the Upscaler which has dune hard-coded as
its solver.
Thank you for introducing yourself and your project, it is always good to know
what people are doing with OPM!
Which brings me to some background on the issue I want to discuss. Support for
more linear solver libraries means we need some control of the data
representation. The upscaler uses Dune's data structures both for building the
problem as well as solving it, but these to me seems unsuited as they're
Dune-specific and hard to convert to some more suited representation.
This is true. I have considered replacing the direct use of Dune::BCRSMatrix in
the upscaler's flow solver (class IncompFlowSolverHybrid) with using
opm/core/linalg/LinearSolverInterface, I have not done so, mostly for the
following reasons:
- Why change working code?
- No resources (time) to do it.
Cue the poor master's student.
And we appreciate it!
- Potential downsides (overhead from matrix conversion etc.)
True, but this is a cost we get by using multiple solvers anyways. It is still
somewhat cheap, and memory bound in most cases.
Yes, so I look forward to seeing what you can manage with this.
[…]
The last week I've been working on a C++ SparseMatrix implementation. It's
written in pure (templated) C++. To make the interface nicer I've also split
the matrix building process into its own class, which solves the state mess
that is prevalent both in struct CSRMatrix and the Dune *CSR matrix classes.
Under the hood it is implemented as plain C arrays, but with the RAII idiom in
place. You can see the WIP in my csr-matrix branch
https://github.com/jorgekva/opm-core/tree/csr-matrix with an emphasis on WIP.
I've constructed the class so that other backends easily can be written and
used, both in construction (e.g. list of lists, dictionaries etc) as well as
representation (BCSR comes to mind).
What I propose is that the call_*.c files are removed in favour of C++
implementations. I also propose that my SparseMatrix gets a place alongside
SparseVector and SparseTable in core/utiliy so we get a uniform and flexible
CSR implementation we can program for that is easily convertible to whatever
format the linear solvers need.
I do not think that the call_*.c files need to be removed, but I do not think
new solver interfaces should be required to provide a C interface.
I don't follow. An extra not, maybe?
The "but" should be a "however".
I would certainly consider such a proposed SparseMatrix class. It should be
placed in opm/core/linalg I think. However, I am afraid that it
may end up as "Yet Another Sparse Matrix" and only add to any complexity or
confusion there already is.
This is what I'm afraid of. linalg vs utility was a rather arbitrary choice -
linalg is fine, I just went with utility because that's where
Sparse{Vector,Table} already was.
Abandoning the Eigen SparseMatrix for opm-autodiff will require a substantial
effort, that I do not think we'll make (at this point at least), so I think
that will remain.
I consider this an issue for the future anyways.
To eventually replace Eigen, we would need to essentially reimplement
large parts of it (not just construction, but all kinds of operators and
utilities), a better choice may then perhaps be to use Eigen as our
adopted representation.
If your proposal gives us greater flexibility for the upscaling code without
hurting performance, I am inclined to accept it anyway, even if that
is the only use case for the new class. As a side note, keep in mind that for
the moment, we require all code to build with g++ 4.4, which places
some limits on C++11 features you can use (for example auto is ok, but lambdas
are not).
I try to limit C++11 features to the ones supported in 4.4, yes. No lambdas, so
I mostly stick with unique_ptr, move semantics, ranged-for and auto.
Ranged for is not in 4.4, the rest should be fine.
Atgeirr
_______________________________________________
Opm mailing list
[email protected]
http://www.opm-project.org/mailman/listinfo/opm
