September 2021 - petsc-users - mailman3.cels.anl.gov

PETSC installation on Cray
by Enrico 13 Sep '21

13 Sep '21

Hi, I'm having some problems installing PETSC with Cray compiler. I use this configuration: ./configure --with-cc=cc --with-cxx=CC --with-fc=0 --with-debugging=1 --with-shared-libraries=1 COPTFLAGS=-O0 CXXOPTFLAGS=-O0 and when I do make all I get the following error because of cmathcalls.h: CC-1043 craycc: ERROR File = /usr/include/bits/cmathcalls.h, Line = 55 _Complex can only be used with floating-point types. __MATHCALL (cacos, (_Mdouble_complex_ __z)); ^ Am I doing something wrong? Regards, Enrico Degregori

4 5

p4est error at NERSC
by Mark Adams 12 Sep '21

12 Sep '21

2 1

DMDA matrices with one sided stencils
by Alfredo J Duarte Gomez 11 Sep '21

11 Sep '21

Good afternoon, I have developed and validated some matrix operators using petsc with a structured dmda. Some of these operators use one-sided stencils at the boundaries, which following the way the dmda uses the stencil width value, requires me to increase the stencil width to accommodate more entries at the boundary only if I want to avoid errors with default options. This is very wasteful and affects my performance, since there are a lot of extra zeros corresponding to the inner points. What is the best way to improve this? I have read in some public threads the possibility of using MatOption to allow us to put more entries into the matrix, but that does not allow me to use MatSetStencil? Alternatively, is there any way to use a larger stencil width and then trim the zero entries that were entered automatically? If there are any other solutions for this problem, please let me know. Thank you, -Alfredo Duarte -- Alfredo Duarte Graduate Research Assistant The University of Texas at Austin

3 2

Combine DM and FFTW
by Gerhard Ungersbäck 11 Sep '21

11 Sep '21

Hello, I want to solve a time dependent differential equation in 3D (Scalar field theory in Hamilton formulation) . The crucial part is that at some time steps I need to FFT the 3D grid. I have written a sequential code without petsc and now I would like to use petsc to get a parallel version. I worked through the examples and now I understand DMDACreate and also the FFTW examples. What is missing though is how I combine both! DMDACreate takes care of ghost points and periodic boundary conditions. As far as I know FFTW requires each process to have a slab of the grid to work. I know how to create this grid with DMDACreate. Normally I would proceed by creating a global vector by DMCreateGlobalVector. But this vector needs then to be linked with fftw arrays. How does this work? Or should I first allocate local memory for fftw and then somehow stitch them together to form a global petsc vector? Thanks best, gerhard

2 1

Petsc with Visual Studio
by Hamid 10 Sep '21

10 Sep '21

3 2

Re: [petsc-users] Ainsworth formula to solve saddle point problems / preconditioner for shell matrices
by Barry Smith 10 Sep '21

10 Sep '21

Olivier, I believe the code is ready for your use. MatInvertVariableBlockEnvelope() will take an MPIAIJ matrix (this would be your C*C' matrix, determine the block diagonal envelope automatically and construct a new MPIAIJ matrix that is the inverse of the block diagonal envelope (by inverting each of the blocks). For small blocks, as is your case, the inverse should inexpensive. I have a test code src/mat/tests/ex178.c that demonstrates its use and tests if for some random blocking across multiple processes. Please let me know if you have any difficulties. The branch is the same name as before barry/2020-10-08/invert-block-diagonal-aij but it has been rebased so you will need to delete your local copy of the branch and then recreate it after git fetch. Barry > On Sep 8, 2021, at 3:09 PM, Olivier Jamond <olivier.jamond(a)cea.fr> wrote: > > Ok thanks a lot! I look forward to hearing from you! > > Best regards, > Olivier > > On 08/09/2021 20:56, Barry Smith wrote: >> >> Olivier, >> >> I'll refresh my memory on this and see if I can finish it up. >> >> Barry >> >> >>> On Sep 2, 2021, at 12:38 PM, Olivier Jamond <olivier.jamond(a)cea.fr <mailto:[email protected]>> wrote: >>> >>> Dear Barry, >>> >>> First I hope that you and your relatives are doing well in these troubled times... >>> >>> I allow myself to come back to you about the subject of being able to compute something like '(C*Ct)^(-1)*C', where 'C' is a 'MPC' matrix that is used to impose some boundary conditions for a structural finite element problem: >>> >>> [K C^t][U]=[F] >>> [C 0 ][L] [D] >>> >>> as we discussed some time ago, I would like to solve such a problem using the Ainsworth method, which involves this '(C*Ct)^(-1)*C'. >>> >>> You kindly started some developments to help me on that, which worked as a 'proof of concept' in sequential, but not yet in parallel, and also kindly suggested that you could extend it to the parallel case (MR: https://gitlab.com/petsc/petsc/-/merge_requests/3544 <https://gitlab.com/petsc/petsc/-/merge_requests/3544>). Can this be still 'scheduled' on your side? >>> >>> Sorry again to "harass" you about that... >>> >>> Best regards, >>> Olivier Jamond >>> >>> On 03/02/2021 08:45, Olivier Jamond wrote: >>>> Dear Barry, >>>> >>>> I come back to you about this topic. As I wrote last time, this is not a "highly urgent" subject (whereas we will have to deal with it in the next months), but it is an important one (especially since the code I am working on just raised significantly its ambitions). So I just would like to check with you that this is still 'scheduled' on your side. >>>> >>>> I am sorry, I feel a little embarrassed about asking you again about your work schedule, but I need some kind of 'visibility' about this topic which will become critical for our application. >>>> >>>> Many thanks for helping me on that! >>>> Olivier >>>> >>>> On 02/12/2020 21:34, Barry Smith wrote: >>>>> >>>>> Sorry I have not gotten back to you quicker, give me a few days to see how viable it is. >>>>> >>>>> Barry >>>>> >>>>> >>>>>> On Nov 25, 2020, at 11:57 AM, Olivier Jamond <olivier.jamond(a)cea.fr <mailto:[email protected]>> wrote: >>>>>> >>>>>> Hi Barry, >>>>>> >>>>>> I come back to you about the feature to unlock the Ainsworth method for saddle point problems in parallel. If I may ask (sorry for that...), is it still on your schedule (I just checked the branch, and it seems 'stuck')? >>>>>> >>>>>> This is not "highly urgent" on my side, but the ability to handle efficiently saddle point problems with iterative solvers will be a critical point for the software I am working on... >>>>>> >>>>>> Many thanks (and sorry again for asking about your work schedule...)! >>>>>> Olivier >>>>>> >>>>>> On 12/10/2020 16:49, Barry Smith wrote: >>>>>>> >>>>>>> >>>>>>>> On Oct 12, 2020, at 6:10 AM, Olivier Jamond <olivier.jamond(a)cea.fr <mailto:[email protected]>> wrote: >>>>>>>> >>>>>>>> Hi Barry, >>>>>>>> >>>>>>>> Thanks for this work! I tried this branch with my code and sequential matrices on a small case: it does work! >>>>>>>> >>>>>>>> >>>>>>> Excellant. I will extend it to the parallel case and get it into our master release. >>>>>>> >>>>>>> We'd be interested in hearing about your convergence and timing experiences when you run largish jobs (even sequentially) since this type of problem comes up relatively frequently and we do need a variety of solvers that can handle it while currently we do not have great tools for it. >>>>>>> >>>>>>> Barry >>>>>>> >>>>>>>> Thanks a lot, >>>>>>>> Olivier >>>>>>>> >>>>>>>> On 09/10/2020 03:50, Barry Smith wrote: >>>>>>>>> >>>>>>>>> Olivier, >>>>>>>>> >>>>>>>>> The branch barry/2020-10-08/invert-block-diagonal-aij contains an example src/mat/tests/ex178.c that shows how to compute inv(CC'). It works for SeqAIJ matrices. >>>>>>>>> >>>>>>>>> Please let us know if it works for you and then I will implement the parallel version. >>>>>>>>> >>>>>>>>> Barry >>>>>>>>> >>>>>>>>> >>>>>>>>>> On Oct 8, 2020, at 3:59 PM, Barry Smith <bsmith(a)petsc.dev <mailto:[email protected]>> wrote: >>>>>>>>>> >>>>>>>>>> >>>>>>>>>> Olivier >>>>>>>>>> >>>>>>>>>> I am working on extending the routines now and hopefully push a branch you can try fairly soon. >>>>>>>>>> >>>>>>>>>> Barry >>>>>>>>>> >>>>>>>>>> >>>>>>>>>>> On Oct 8, 2020, at 3:07 PM, Jed Brown <jed(a)jedbrown.org <mailto:[email protected]>> wrote: >>>>>>>>>>> >>>>>>>>>>> Olivier Jamond <olivier.jamond(a)cea.fr <mailto:[email protected]>> writes: >>>>>>>>>>> >>>>>>>>>>>>> Given the structure of C it seems you should just explicitly construct Sp and use GAMG (or other preconditioners, even a direct solver) directly on Sp. Trying to avoid explicitly forming Sp will give you a much slower performing solving for what benefit? If C was just some generic monster than forming Sp might be unrealistic but in your case CCt is is block diagonal with tiny blocks which means (C*Ct)^(-1) is block diagonal with tiny blocks (the blocks are the inverses of the blocks of (C*Ct)). >>>>>>>>>>>>> >>>>>>>>>>>>> Sp = Ct*C + Qt * S * Q = Ct*C + [I - Ct * (C*Ct)^(-1)*C] S [I - Ct * (C*Ct)^(-1)*C] >>>>>>>>>>>>> >>>>>>>>>>>>> [Ct * (C*Ct)^(-1)*C] will again be block diagonal with slightly larger blocks. >>>>>>>>>>>>> >>>>>>>>>>>>> You can do D = (C*Ct) with MatMatMult() then write custom code that zips through the diagonal blocks of D inverting all of them to get iD then use MatPtAP applied to C and iD to get Ct * (C*Ct)^(-1)*C then MatShift() to include the I then MatPtAP or MatRAR to get [I - Ct * (C*Ct)^(-1)*C] S [I - Ct * (C*Ct)^(-1)*C] then finally MatAXPY() to get Sp. The complexity of each of the Mat operations is very low because of the absurdly simple structure of C and its descendants. You might even be able to just use MUMPS to give you the explicit inv(C*Ct) without writing custom code to get iD. >>>>>>>>>>>> >>>>>>>>>>>> At this time, I didn't manage to compute iD=inv(C*Ct) without using >>>>>>>>>>>> dense matrices, what may be a shame because all matrices are sparse . Is >>>>>>>>>>>> it possible? >>>>>>>>>>>> >>>>>>>>>>>> And I get no idea of how to write code to manually zip through the >>>>>>>>>>>> diagonal blocks of D to invert them... >>>>>>>>>>> >>>>>>>>>>> You could use MatInvertVariableBlockDiagonal(), which should perhaps return a Mat instead of a raw array. >>>>>>>>>>> >>>>>>>>>>> If you have constant block sizes, MatInvertBlockDiagonalMat will return a Mat. >>>>>>>>>> >>>>>>>>> >>>>>>> >>>>> >>

1 0

this is what mangers are for ...
by Mark Adams 10 Sep '21

10 Sep '21

Hi Sherry, Here is a request that you do not see every day. Can you please verify that the back of your badge says this for me? I need a badge for access to Fugaku and my badge is old and unreadable. They apparently grabbed this text from some other LBNL badge and want to verify. Thanks, Mark "EMERGENCY STATUS ANNOUNCEMENT IF FOUND,DROP IN ANY MAILBOX,POSTMASTER,POSTAGE GUARANTEED. RETURN TO: LAWRENCE BERKELEY NATIONAL LABORATORY ONE CYCLOTRON ROAD BERKELEY, CALIFORNIA 94720 THIS CREDENTIAL IS THE PROPERTY OF THE U.S. GOVERNMENT. THE COUNTERFEIT, ALTERATION, OR MISUSE IS A VIOLATION OF SECTION 499 AND 701, TITLE 18,UNITED STATES CODE. OPERATED UNDER CONTRACT NO. DE-AC03-76SF00098 WITH THE U.S.DEPARTMENT OF ENERGY"

1 1

petsc_with_c++_python
by govind sharma 08 Sep '21

08 Sep '21

Hi, I need a 2d poisson solver such that I can use petsc objects between C++ and python interactively. I explain this way: -->1 Let's 2D poisson solver written in C++ with petsc -->2 With python notebook interact with solve in parallel with petscpy --> IPython parallel I have seen an example of doing this type by Aron group with C but I find it incomprehensible to understand. Can I get some clarification on this? How do I proceed? Regards, Govind

2 3

Mat preallocation in case of variable stencils
by Matteo Semplice 07 Sep '21

07 Sep '21

Hi. We are writing a code for a FD scheme on an irregular domain and thus the local stencil is quite variable: we have inner nodes, boundary nodes and inactive nodes, each with their own stencil type and offset with respect to the grid node. We currently create a matrix with DMCreateMatrix on a DMDA and for now have set the option MAT_NEW_NONZERO_LOCATIONS to PETSC_TRUE, but its time to render the code memory-efficient. The layout created automatically is correct for inner nodes, wrong for boundary ones (off-centered stencils) and redundant for outer nodes. After the preprocessing stage (including stencil creation) we'd be in position to set the nonzero pattern properly. Do we need to start from a Mat created by CreateMatrix? Or is it ok to call DMCreateMatrix (so that the splitting among CPUs and the block size are set by PETSc) and then call a MatSetPreallocation routine? Also, I've seen in some examples that you call the Seq and the MPI preallocation routines in a row. Does this work because MatMPIAIJSetPreallocation silently does nothing on a Seq matrix and viceversa? Thanks Matteo

4 5

Slow convergence while parallel computations.
by Наздрачёв Виктор 06 Sep '21

06 Sep '21

Dear all, I have a 3D elasticity problem with heterogeneous properties. There is unstructured grid with aspect ratio varied from 4 to 25. Zero Dirichlet BCs are imposed on bottom face of mesh. Also, Neumann (traction) BCs are imposed on side faces. Gravity load is also accounted for. The grid I use consists of 500k cells (which is approximately 1.6M of DOFs). The best performance and memory usage for single MPI process was obtained with HPDDM(BFBCG) solver and bjacobian + ICC (1) in subdomains as preconditioner, it took 1 m 45 s and RAM 5.0 GB. Parallel computation with 4 MPI processes took 2 m 46 s when using 5.6 GB of RAM. This because of number of iterations required to achieve the same tolerance is significantly increased. I`ve also tried PCGAMG (agg) preconditioner with ICС (1) sub-precondtioner. For single MPI process, the calculation took 10 min and 3.4 GB of RAM. To improve the convergence rate, the nullspace was attached using MatNullSpaceCreateRigidBody and MatSetNearNullSpace subroutines. This has reduced calculation time to 3 m 58 s when using 4.3 GB of RAM. Also, there is peak memory usage with 14.1 GB, which appears just before the start of the iterations. Parallel computation with 4 MPI processes took 2 m 53 s when using 8.4 GB of RAM. In that case the peak memory usage is about 22 GB. Are there ways to avoid decreasing of the convergence rate for bjacobi precondtioner in parallel mode? Does it make sense to use hierarchical or nested krylov methods with a local gmres solver (sub_pc_type gmres) and some sub-precondtioner (for example, sub_pc_type bjacobi)? Is this peak memory usage expected for gamg preconditioner? is there any way to reduce it? What advice would you give to improve the convergence rate with multiple MPI processes, but keep memory consumption reasonable? Kind regards, Viktor Nazdrachev R&D senior researcher Geosteering Technologies LLC

7 15