<div dir="ltr"><div><div>Hey Oli,<br><br></div>as I was implementing the new matrix build mode back then, I can give you some numbers and an idea of why it is so important.<br>The way it was implemented before, the following steps were done strictly sequentially<br></div><div>1) Determine sparsity pattern and store it in a suitable (map-like) data structure<br></div><div>2) Allocate storage for the BCRSMatrix<br></div><div>3) Build the sparsity pattern on the matrix using the data from 1)<br><br></div><div>The key problem here is peek memory usage: Step 1) has memory complexity n*log(n) (n being the number of matrix entries),<br>Step 2) has complexity n. Both data structures need to be kept around until step 3 is finished. Only then, the data structure<br>from 1) can be discarded.<br><br></div><div>The new algorithm works more like this:<br></div><div>1) Get some pattern statistics from the user<br></div><div>2) Allocate storage for the BCRSMatrix based on the given statistics<br></div><div>3) Build the sparsity pattern on-the-fly directly into the allocated storage<br></div><div>4) "Compress" the matrix (remove holes from underfill etc.)<br></div><div><br></div><div>Given a vanilla 2D Q1 poisson problem with around a million unknowns, the peak usage dropped from 3.6 GB to 1.1 GB.<br></div><div>You can definitely see that the log(n) factor on the pattern-map's complexity really does contribute here. You might<br></div><div>also extrapolate that for even larger problems (both in number of dofs and in entries per row), the problem potentially poses a<br></div><div>bottleneck for the entire simulation code.<br><br></div><div>I agree, that some good documentation on the why and how of this would be important.<br></div><div><br></div><div>Best,<br></div><div>Dominic<br></div><div><br><br></div></div><div class="gmail_extra"><br><div class="gmail_quote">On Mon, Jul 20, 2015 at 2:50 PM, Steffen Müthing <span dir="ltr"><<a href="mailto:steffen.muething@iwr.uni-heidelberg.de" target="_blank">steffen.muething@iwr.uni-heidelberg.de</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hi Oliver,<br>
<span class=""><br>
> Am 13.07.2015 um 11:58 schrieb Oliver Sander <<a href="mailto:sander@igpm.rwth-aachen.de">sander@igpm.rwth-aachen.de</a>>:<br>
><br>
> Dear pdelab,<br>
><br>
> I just switched from the old ISTLMatrixBackend to the new istl::BCRSMatrixBackend.<br>
> Everything appears to work nicely, but I find the new implementation more cumbersome<br>
> to use.<br>
><br>
> Before, I could simple instantiate GridOperator with the correct matrix<br>
> backend type, and the corresponding object was created automatically for me.<br>
> With the new interface, I have to create a backend object myself (that's acceptable),<br>
> and additionally I have to give it information about the expected matrix density.<br>
> But I am solving small problems on small grids. Memory consumption and wall-clock<br>
> time are not important to me. Therefore I don't want to be forced to think about<br>
> matrix occupation patterns. Is it possible to make the backend object default<br>
> constructible again (maybe choosing 3^d as a default value), so people like me<br>
> don't have to think about issues that only affect HPC projects?<br>
<br>
</span>The new construction method is important not just for HPC, but also for people just using PDELab<br>
on a small scale (laptops etc.). I don’t have the exact numbers right now (Dominic did the benchmarks<br>
back when this was implemented), but for a simple Poisson problem with P1 elements, the new method<br>
reduces the maximum memory consumption by a factor of 2 - 3 and can cut down the run time (of the<br>
overall program) by an order of magnitude for larger problems. Unfortunately, it’s really not possible to<br>
specify a sensible global default (your idea of using 3^d is e.g. way too small for non-blocked DG spaces,<br>
where a good formula would be (average_number_of_neighbor_intersections + 1) * (k+1)^d ).<br>
<br>
That said, we really need to document this a lot better somewhere and give people some good "cooking<br>
recipes" for calculating these numbers. I’ll open a bug for that…<br>
<br>
Steffen<br>
<br>
><br>
> Thanks,<br>
> Oliver<br>
><br>
><br>
> _______________________________________________<br>
> dune-pdelab mailing list<br>
> <a href="mailto:dune-pdelab@dune-project.org">dune-pdelab@dune-project.org</a><br>
> <a href="http://lists.dune-project.org/mailman/listinfo/dune-pdelab" rel="noreferrer" target="_blank">http://lists.dune-project.org/mailman/listinfo/dune-pdelab</a><br>
<br>
<br>_______________________________________________<br>
dune-pdelab mailing list<br>
<a href="mailto:dune-pdelab@dune-project.org">dune-pdelab@dune-project.org</a><br>
<a href="http://lists.dune-project.org/mailman/listinfo/dune-pdelab" rel="noreferrer" target="_blank">http://lists.dune-project.org/mailman/listinfo/dune-pdelab</a><br>
<br></blockquote></div><br></div>