[Dune] Preconditioner

Mike Rosing eresrch at sdf.org
Tue Sep 9 17:00:29 CEST 2025


Thanks!  I will go take a look at both dune-common and dune-functions.
The switchCases looks interesting, and will be good to help me learn
more C++.

Mike

On Tue, 9 Sep 2025, Simon Praetorius wrote:

> Date: Tue, 9 Sep 2025 15:24:13 +0200
> From: Simon Praetorius <simon.praetorius at tu-dresden.de>
> To: dune at lists.dune-project.org
> Subject: Re: [Dune] Preconditioner
> 
> Hi,
>
> Perhaps I can contribute a few more thoughts: if the size of the blocks is 
> variable but fixed, in the sense that each setup requires a different size, 
> but each block is the same size and not changed within a setup, then 
> BCRSMatrix<FieldMatrix> would be the best idea, also in terms of performance. 
> You could implement a simple switch-case branch in the main program to define 
> the block size. Then you can map everything in one program. There is also a 
> helper structure for this in dune-common:
>
> Dune::Hybrid::switchCases(std::index_sequence<1,2,3,4,5,...>{}, n, [&](auto 
> N) { run<N>(); });
>
> Here, all possible sizes n in {1,2,3,4,5...} are tested and the correct one 
> is then passed as a static constant to the passed (lambda) function, in which 
> you can then use a FieldMatrix<double,N,N>.
>
> An alternative you mentioned is to flatten the matrix structure. There are 
> various ways to do this. I use dune-functions, where blocking can be defined 
> in the function space basis using a so-called index merging strategy. This 
> allows you to map an index of a local basis function on a grid element to a 
> global multi-index or flat index, which you can then use to address the 
> matrix entries. With the help of the MatrixIndexSet, you can use this to 
> build the non-zero structure of the BCRSMatrix accordingly. See the examples 
> in dune-functions for an explanation how to do this.
>
> Best,
> Simon
>
>
> Am 09.09.25 um 14:59 schrieb Mike Rosing:
>> Howdy Christian,
>> 
>> Recompiling everytime I change size is annoying, but it's also not
>> that big a deal. Worst case if I want different sizes I can just use
>> different named variables - so I can see how to get there from here.
>> 
>> I just have to change my thought process to the reality of the code.
>> 
>> I did try flattening the matrix to BCRSMatrix<double> but I must have
>> done something wrong because DUNE aborted with a NaN error. So I think
>> being safe is the way to go.
>> 
>> Thanks!
>> Mike
>> 
>> On Tue, 9 Sep 2025, Christian Engwer wrote:
>> 
>>> Date: Tue, 09 Sep 2025 07:55:16 +0200
>>> From: Christian Engwer <christian.engwer at uni-muenster.de>
>>> To: dune at lists.dune-project.org, Mike Rosing <eresrch at sdf.org>
>>> Subject: Re: [Dune] Preconditioner
>>> 
>>> Dear Mike,
>>> 
>>> I'm afraid also we DUNE people need to dig deeper into this.
>>> 
>>> The issue here is that you're using Matrix for the blocks, which have 
>>> dynamic size. The usual pattern is to use FieldMatrix, which has static 
>>> size and a default constructor. The algorithms use this to decide where to 
>>> invest directly.
>>> 
>>> The genetic algorithms usually expect BCRSMatrix to contain dense data, 
>>> but Matrix dies not where that everything below is dense.
>>> 
>>> The safe way would be to use BCRSMatrix<FieldMatrix<double,n,n>> and 
>>> BlockVector<FieldVector<double, n>>
>>> 
>>> Do you really need the dynamic size?
>>> 
>>> I didn't test myself, but perhaps BCRSMatrix<DynamicMatrix<double,n,n>> 
>>> and BlockVector<DynamicVector<double, n>> work. I would expect relaxation 
>>> methods to work, but not AMG. Für ILU it might be a problem, that the 
>>> constructor of DynamicMatrix needs the size information.
>>> 
>>> Best
>>> Christian
>>> 
>>> Am 6. September 2025 21:04:48 MESZ schrieb Mike Rosing <eresrch at sdf.org>:
>>>> I figured out how to get the MatrixAdaptor to be happy with my 
>>>> structures:
>>>> 
>>>> using Occupy = Matrix<double>;
>>>> using ElementVec = BlockVector<double>;
>>>> 
>>>>  MatrixAdapter<BCRSMatrix<Occupy>, BlockVector<ElementVec>,
>>>>                    BlockVector<ElementVec>> 
>>>> linearOperator(stiffns);
>>>> 
>>>> but none of the preconditioners are happy with this.
>>>> The ILU preconditioner gives the attached error output.
>>>> 
>>>> Does it make any sense to attempt to create a preconditioner that can 
>>>> deal with this structure, or should I just copy all the data into a flat 
>>>> matrix and vector with no sub-structure? Because I haven't gotten to call 
>>>> the solver yet, so this might be a lot deeper than I want to know about.
>>>> 
>>>> Thanks,
>>>> 
>>>> 
>>>> Mike Rosing
>>>> Engineering Research
>>>> Madison WI USA
>>>> Elliptic Curve Crypto:  http://mng.bz/D9NA
>>> 
>> 
>> Mike Rosing
>> Engineering Research
>> Madison WI USA
>> Elliptic Curve Crypto:  http://mng.bz/D9NA
>> 
>> _______________________________________________
>> Dune mailing list
>> Dune at lists.dune-project.org
>> https://lists.dune-project.org/mailman/listinfo/dune
>
> -- 
> Dr. Simon Praetorius
> Technische Universität Dresden
> Institute of Scientific Computing
> phone: +49 351 463-34432
> mail: simon.praetorius@​tu-dresden.de
>
>
> _______________________________________________
> Dune mailing list
> Dune at lists.dune-project.org
> https://lists.dune-project.org/mailman/listinfo/dune
>
>

Mike Rosing
Engineering Research
Madison WI USA
Elliptic Curve Crypto:  http://mng.bz/D9NA


More information about the Dune mailing list