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

Sum of the absolute values of each row's components
by Salazar De Troya, Miguel 24 May '21

24 May '21

Hello, I am simply interested in obtaining a vector with x_i = \sum_j | A_{i, j} | for each row “i” in the matrix A. I found MatGetColumnNorms(), but no row version. I am wondering if it is more efficient to calculate the transpose A^T and then MatGetColumnNorms() or maybe iterate through each row with MatGetRow() and calculate \sum_j | A_{i, j} | by hand by myself. Thanks Miguel Miguel A. Salazar de Troya Postdoctoral Researcher, Lawrence Livermore National Laboratory B141 Rm: 1085-5 Ph: 1(925) 422-6411

2 2

MatChop
by Emmanuel Ayala 21 May '21

21 May '21

Hi everybody, I just updated petsc from version 13 to 15. Before the update everything works well, then my code give me an error: [9]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [9]PETSC ERROR: Invalid argument [9]PETSC ERROR: Setting off process row 53484 even though MatSetOption(,MAT_NO_OFF_PROC_ENTRIES,PETSC_TRUE) was set [9]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting. [9]PETSC ERROR: Petsc Release Version 3.15.0, Mar 30, 2021 [9]PETSC ERROR: ./comp on a arch-linux-c-opt-O2-superlud_mumps_hyp named eayala by ayala Fri May 21 10:40:36 2021 [9]PETSC ERROR: Configure options --with-debugging=0 COPTFLAGS="-O2 -march=native -mtune=native" CXXOPTFLAGS="-O2 -march=native -mtune=native" FOPTFLAGS="-O2 -march=native -mtune=native" --download-mpich --download-hypre --download-mumps --download-scalapack --download-parmetis --download-metis --download-superlu_dist --download-cmake --download-fblaslapack=1 --with-cxx-dialect=C++11 The error appears after a matrix assembly, the matrix was created with DMCreateMatrix and updated with MatSetValuesLocal. A Little work around I found the solution, avoid using MatChop on this matrix, but I still need to use MatChop. There is any reason to have this problem? Thanks in advance.

5 9

Nested SNES in FormFunction
by Pierre Bernigaud 21 May '21

21 May '21

Greetings, I am currently working on a CFD solver using PETSc. I have a non linear system which is solved using 2D_DMDA/SNES, and submitted to boundary conditions that are treated implicitly and updated in the FormFunction. The calculation of one of these boundary conditions requires the resolution of an other non linear system. I am hence using a nested 1D_DMDA/SNES system within the FormFunction of my main SNES solver to solve for this boundary condition. This is working fine, but doing a scalability study we found out that this causes the code to show sub-par acceleration properties. Have you ever encountered this kind of nested SNES application, and are there some critical points to be aware of in order to avoid a loss of performances? For instance, the sub 1D_DMDA/SNES objects are created and destroyed at each update of the boundary, hence at each call to FormFunction, which results in an important number of object manipulation. Could this be a problem? Furthermore, the use of a sub 1D_DMDA/SNES allows to use multiple processors to solve for the boundary condition, composed of say N cells. When running the code with M > N processors, everything is working great, but I am curious about the state of the (M-N) processors which aren’t working on boundary condition problem. Do they just stay idle? Thank you for your help. Respectfully, Pierre Bernigaud

2 1

nghosts in 'DMMoabLoadFromFile'
by Bhamidipati, Vikram 20 May '21

20 May '21

Hello, I am in the process of updating PetSc versions and I see that 'DMMoabLoadFromFile' function (in dmmutil.cxx) has a new argument 'nghosts'. For those of us who don't use ghost cells should we set it to 0? Thanks, Vikram --------------------------------------------------- Vikram Bhamidipati Senior Research Engineer Computational Material Integrity Section Materials Engineering Department Mechanical Engineering Division Ph: (210) 522-2576

3 3

DMPlex and Boundary facets
by Karin&NiKo 19 May '21

19 May '21

Dear PETSc team, I have tried to load a test mesh available in Gmsh' s demos directory (share/doc/gmsh/demos/simple_geo/filter.geo, attached to this email) as a DMPlex. So I produced a msh4 file by doing : gmsh -3 filter.geo -o /tmp/test.msh4 Then I used src/dm/impls/plex/tutorials/ex2.c to load the mesh by doing : ./ex2 -filename /tmp/test.msh4 Unfortunately I get the error : [0]PETSC ERROR: --------------------- Error Message -------------------------------------------------------------- [0]PETSC ERROR: No support for this operation for this object type [0]PETSC ERROR: Could not determine Plex facet for Gmsh element 1268 (Plex cell 12681) The error seems to come from the fact that the msh file contains tets *and* facets *only on the Physical entities* (aka parts of the mesh boundary where the user will assign Dirichlet or Neuman conditions). If I suppress these facets by commenting the "Physical Surface" lines in the geo file and regenerating the mesh, everything is fine. But the use of these "Physical" stuff is very common in lots of finite element codes in order to assign boundary conditions. How should I do to keep these boundary groups of 2D elements (with corresponding names) ? Thanks for your help, Nicolas

3 11

Using SNES to solve Non linear elasticity, passing matrix to form function
by Kaushik Vijaykumar 19 May '21

19 May '21

Hi everyone, We are in the process of integrating nonlinear solution to our FE code through SNES. The first steps that I understood that need to be done is to be able to pass the assembled stiffness matrix K and force vector F to the " *formfunction*" to calculate the residual and "*formjacobian*" to calculate the tangent stiffness matrix. To do so, I have defined a struct: *extern PetscErrorCode FormJacobian(SNES,Vec,Mat,Mat,void*);extern PetscErrorCode FormFunction(SNES,Vec,Vec,void*);* // define a struct type to pass K and f as "user context" *typedef struct { Mat K; Vec f; } K_and_f;* In the main program, the struct is declared *int main()* *{* *K_and_f main_K_and_f;* // declare the struct // SNES - Populate K and f into the struct *main_K_and_f.K = K;* // K matrix *main_K_and_f.f = f; *// f vector .... .... *}* In form function *PetscErrorCode FormFunction(SNES snes, Vec x, Vec F, void* ctx) { PetscErrorCode ierr; PetscReal *ax,*c; PetscReal *aF; PetscScalar *Kc,v; PetscInt nlocal,m,n,i,j,index[100000]; * // Create local Vec and Mat * Mat Kloc;* * Vec Floc, Uloc, KUloc, resloc;* * K_and_f* ptr = (K_and_f*) ctx; *// cast the pointer to void into pointer to struct // Get local F array, FLOC * ierr = VecGhostGetLocalForm(ptr->f,&Floc);CHKERRQ(ierr);* I am able to get the f array from the main program using " VecGhostGetLocalForm" and the vec is correct. However, I am having trouble with reading the K matrix in formfunction. The following trial gives me incorrect K values in Kloc: *ierr = MatSeqAIJGetArray(ptr->K,&Kc);CHKERRQ(ierr); * *ierr = MatGetLocalSize(ptr->K,&m,&n); CHKERRQ(ierr);* *for(i=0;i<m;i++) { for(j=0;j<n;j++) { v = *(Kc + (i*n+j)); ierr = MatSetValues(Kloc,1,&i,1,&j,&v,INSERT_VALUES);CHKERRQ(ierr); }* * } * When I compare K and Kloc, they are not identical Please let me know if there is an equivalent function like *VecGhostGetLocalForm()* for Matrices. If not, is there a better way to do this. Any guidance/help is greatly appreciated. Thanks Kaushik

4 7

The method of print the preconditioned matrix
by 杨念 19 May '21

19 May '21

Dear Sir or Madam:   Hello, I am a new PETSc learner. I'm having some problem with it and I'm using version 3.6.3 of PETSc. I can only print out the Jacobi matrix (J) with -mat_view from the command line, but I also need to print out the preconditioning operator (M^(-1)) or the preconditioned matrix (JM^(-1)). Is there a function like that in PETSC? Could you please give me some suggestions? Thank you so much.   Best regards,   Nian Yang

2 1

DMPlex from periodic Gmsh and ghost cells
by Thibault Bridel-Bertomeu 19 May '21

19 May '21

Dear all, I am playing around with creating DMPlex from a periodic Gmsh mesh (still in the same finite volume code that solves the Euler equations) because I want to run some classical periodic test cases from the literature (like convecting a vortex and so on). I got from the mailing list and the examples that I gotta use -dm_plex_gmsh_periodic for the reader and the DM creator to do their job properly and read the $Periodic section from the Gmsh mesh. That works fine, I can write the DM in HDF5 format afterwards and the periodicity is clearly there (i can see in Paraview that opposite sides are "attached" to each other). Now, in the context of the finite volume solver, I wrote my one right-hand side term computing routine, and I rely heavily on routines such as DMPlexGetFaceFields. Those routines, in turn rely on the ghost cells and ... I am having an issue with DMPlexConstructGhostCells. When I use -dm_plex_gmsh_periodic, DMPlexCreateFromFile works fine, DMSetBasicAdjacency(dm, true, false) works fine, DMPlexDistribute works fine, but when I get to DMPlexConstructGhostCells, I get "DM has boundary face 2048 with 2 support cells" I understand the DM begin periodic, even the "boundary" faces have two neighbors because of the topological folding, so ... yeah, the ghost cells cannot really exist. But then how can I proceed ? Because I need that "ghost" label for all the routines I call from the RHS computing routine ... I attach to this email the periodic Gmsh mesh I am playing with, and the routine from my code where I create the DM and so on: subroutine initmesh PetscErrorCode :: ierr DM :: dmDist, dmGhost integer :: overlap PetscViewer :: vtkViewer call PetscLogEventBegin(timer_initmesh, ierr); CHKERRA(ierr) ! Number of neighbours taken into account in MP communications(1 - Order 1; 2 - Order 2) overlap = 1 call PetscPrintf(PETSC_COMM_WORLD, "Initializing mesh...\n", ierr) ; CHKERRA(ierr) ! Force DMPlex to use gmsh marker ! call PetscOptionsSetValue(PETSC_NULL_OPTIONS, "-dm_plex_gmsh_use_marker", "true", ierr); CHKERRA(ierr) ! Read mesh from file name 'meshname' call DMPlexCreateFromFile(PETSC_COMM_WORLD, meshname, PETSC_TRUE, dm, ierr); CHKERRA(ierr) ! Distribute on processors ! Start with connectivity call DMSetBasicAdjacency(dm, PETSC_TRUE, PETSC_FALSE, ierr) ; CHKERRA(ierr) ! Distribute on processors call DMPlexDistribute(dm, overlap, PETSC_NULL_SF, dmDist, ierr) ; CHKERRA(ierr) ! Security check if (dmDist /= PETSC_NULL_DM) then ! Destroy previous dm call DMDestroy(dm, ierr) ; CHKERRA(ierr) ! Replace with dmDist dm = dmDist end if ! Finalize setup of the object call DMSetFromOptions(dm, ierr) ; CHKERRA(ierr) ! Boundary condition with ghost cells call DMPlexConstructGhostCells(dm, PETSC_NULL_CHARACTER, PETSC_NULL_INTEGER, dmGhost, ierr); CHKERRA(ierr) ! Security check if (dmGhost /= PETSC_NULL_DM) then ! Destroy previous dm call DMDestroy(dm, ierr) ; CHKERRA(ierr) ! Replace with dmGhost dm = dmGhost end if if (debug) then ! Show in terminal call PetscPrintf(PETSC_COMM_WORLD, ":: [DEBUG] Visualizing DM in console ::\n", ierr); CHKERRA(ierr) call DMView(dm, PETSC_VIEWER_STDOUT_WORLD, ierr) ; CHKERRA(ierr) ! VTK viewer call PetscViewerCreate(PETSC_COMM_WORLD, vtkViewer, ierr) ; CHKERRA(ierr) call PetscViewerSetType(vtkViewer, PETSCVIEWERHDF5, ierr) ; CHKERRA(ierr) call PetscViewerFileSetMode(vtkViewer, FILE_MODE_WRITE, ierr) ; CHKERRA(ierr) call PetscViewerFileSetName(vtkViewer, "debug_initmesh.h5", ierr) ; CHKERRA(ierr) call DMView(dm, vtkViewer, ierr) ; CHKERRA(ierr) call PetscViewerDestroy(vtkViewer, ierr) ; CHKERRA(ierr) end if call PetscPrintf(PETSC_COMM_WORLD, "Done !\n", ierr) ; CHKERRA(ierr) call PetscLogEventEnd(timer_initmesh, ierr); CHKERRA(ierr) end subroutine initmesh Thank you very much for your help ! Best regards, Thibault Bridel-Bertomeu

2 5

Re: [petsc-users] petsc-users Digest, Vol 149, Issue 47
by Sajid Ali 18 May '21

18 May '21

You could use VecPermulte and MatPermute from the PETSc API to permute the vectors/matrices. However, MatPermute creates a new matrix and even though VecPermute permutes the vector (locally) in-place, it allocates a temporary array and frees the original array. Since you are working on dense matrices and vectors and desire to avoid temporary allocation, you could use BLAS level 1 swap function (used by VecSwap to swap two vectors) which will probably be the most optimized version for the hardware you're using (since it's implemented by platform specific intrinsics and assembly). ---------------------------------------------------------------------- > > Message: 1 > Date: Tue, 18 May 2021 11:44:36 -0500 > From: Barry Smith <bsmith(a)petsc.dev> > To: Roland Richter <roland.richter(a)ntnu.no> > Cc: PETSc <petsc-users(a)mcs.anl.gov> > Subject: Re: [petsc-users] Efficient FFTShift-implementation for > vectors/matrices > Message-ID: <83012843-A364-4CE1-A92B-1768EB11591E(a)petsc.dev> > Content-Type: text/plain; charset="us-ascii" > > > I found a variety of things on the web, below. I don't understand this > but for the even case it seems one simply modifies the input matrix before > the FFT http://www.fftw.org/faq/section3.html#centerorigin > > > https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-c… > < > https://stackoverflow.com/questions/5915125/fftshift-ifftshift-c-c-source-c… > > > > https://www.dsprelated.com/showthread/comp.dsp/20790-1.php < > https://www.dsprelated.com/showthread/comp.dsp/20790-1.php> > > > > > On May 18, 2021, at 9:48 AM, Roland Richter <roland.richter(a)ntnu.no> > wrote: > > > > Dear all, > > > > I tried to implement the function fftshift from numpy (i.e. swap the > > half-spaces of all axis) for row vectors in a matrix by using the > > following code > > > > void fft_shift(Mat &fft_matrix) { > > PetscScalar *mat_ptr; > > MatDenseGetArray (fft_matrix, &mat_ptr); > > PetscInt r_0, r_1; > > MatGetOwnershipRange(fft_matrix, &r_0, &r_1); > > PetscInt local_row_num = r_1 - r_0; > > arma::cx_mat temp_mat(local_row_num, Ntime, arma::fill::zeros); > > for(int i = 0; i < Ntime; ++i) { > > const PetscInt row_shift = i * local_row_num; > > for(int j = 0; j < local_row_num; ++j) { > > const PetscInt cur_pos = j + row_shift; > > if(i < (int)(Ntime / 2)) > > temp_mat(j, i + int(Ntime / 2)) = *(mat_ptr + > cur_pos); > > else > > temp_mat(j, i - int(Ntime / 2)) = *(mat_ptr + > cur_pos); > > } > > } > > for(int i = 0; i < Ntime; ++i) { > > const PetscInt row_shift = i * local_row_num; > > for(int j = 0; j < local_row_num; ++j) { > > const PetscInt cur_pos = j + row_shift; > > *(mat_ptr + cur_pos) = temp_mat(j, i); > > } > > } > > MatDenseRestoreArray (fft_matrix, &mat_ptr); > > } > > > > but I do not like the approach of having a second matrix as temporary > > storage space. Are there more efficient approaches possible using > > PETSc-functions? > > > > Thanks! > > > > Regards, > > > > Roland Richter > > > >

1 0

configure error
by Mark Adams 18 May '21

18 May '21

I can build a non-kokkos PETSc here (Cori/GPU) but I get this error with Kokkos. Any suggestions? Thanks, Mark ============================================================================================= ============================================================================================= Trying to download git://https://bitbucket.org/petsc/pkg-sowing.git for SOWING ============================================================================================= ============================================================================================= Trying to download git://https://bitbucket.org/petsc/pkg-sowing.git for SOWING ============================================================================================= ============================================================================================= Trying to download https://bitbucket.org/petsc/pkg-sowing/get/v1.1.26-p1.tar.gz for SOWING ============================================================================================= ******************************************************************************* UNABLE to CONFIGURE with GIVEN OPTIONS (see configure.log for details): ------------------------------------------------------------------------------- Error during download/extract/detection of SOWING: Unable to clone sowing Could not execute "['git clone https://bitbucket.org/petsc/pkg-sowing.git /global/u2/m/madams/petsc/arch-cori-gpu-opt-kokkos-gcc/externalpackages/git.sowing']": Cloning into '/global/u2/m/madams/petsc/arch-cori-gpu-opt-kokkos-gcc/externalpackages/git.sowing'... fatal: unable to access 'https://bitbucket.org/petsc/pkg-sowing.git/': The requested URL returned error: 503 Unable to download package SOWING from: git:// https://bitbucket.org/petsc/pkg-sowing.git * If URL specified manually - perhaps there is a typo? * If your network is disconnected - please reconnect and rerun ./configure * Or perhaps you have a firewall blocking the download * You can run with --with-packages-download-dir=/adirectory and ./configure will instruct you what packages to download manually * or you can download the above URL manually, to /yourselectedlocation and use the configure option: --download-sowing=/yourselectedlocation Unable to clone sowing Could not execute "['git clone https://bitbucket.org/petsc/pkg-sowing.git /global/u2/m/madams/petsc/arch-cori-gpu-opt-kokkos-gcc/externalpackages/git.sowing']": Cloning into '/global/u2/m/madams/petsc/arch-cori-gpu-opt-kokkos-gcc/externalpackages/git.sowing'... fatal: unable to access 'https://bitbucket.org/petsc/pkg-sowing.git/': The requested URL returned error: 503 Unable to download package SOWING from: git:// https://bitbucket.org/petsc/pkg-sowing.git * If URL specified manually - perhaps there is a typo? * If your network is disconnected - please reconnect and rerun ./configure * Or perhaps you have a firewall blocking the download * You can run with --with-packages-download-dir=/adirectory and ./configure will instruct you what packages to download manually * or you can download the above URL manually, to /yourselectedlocation and use the configure option: --download-sowing=/yourselectedlocation file could not be opened successfully Downloaded package SOWING from: https://bitbucket.org/petsc/pkg-sowing/get/v1.1.26-p1.tar.gz is not a tarball. [or installed python cannot process compressed files] * If you are behind a firewall - please fix your proxy and rerun ./configure For example at LANL you may need to set the environmental variable http_proxy (or HTTP_PROXY?) to http://proxyout.lanl.gov * You can run with --with-packages-download-dir=/adirectory and ./configure will instruct you what packages to download manually * or you can download the above URL manually, to /yourselectedlocation/v1.1.26-p1.tar.gz and use the configure option: --download-sowing=/yourselectedlocation/v1.1.26-p1.tar.gz *******************************************************************************

3 17