Merge pull request #17 from firedrakeproject/connorjward/merge-upstream

Merge upstream

config/BuildSystem/config/packages/MUMPS.py

@@ -3,7 +3,7 @@
 class Configure(config.package.Package):
   def __init__(self, framework):
     config.package.Package.__init__(self, framework)
-    self.version          = '5.6.1'
+    self.version          = '5.6.2'
     self.minversion       = '5.2.1'
     self.versionname      = 'MUMPS_VERSION'
     self.requiresversion  = 1
@@ -22,6 +22,7 @@ def setupHelp(self, help):
     import nargs
     config.package.Package.setupHelp(self, help)
     help.addArgument('MUMPS', '-with-mumps-serial', nargs.ArgBool(None, 0, 'Use serial build of MUMPS'))
+    help.addArgument('MUMPS', '-download-mumps-avoid-mpi-in-place', nargs.ArgBool(None, 0, 'Let MUMPS not use MPI_IN_PLACE. Since MUMPS-5.6.2, it can be used to avoid a bug in MPICH older than 4.0b1'))
     return
 
   def setupDependencies(self, framework):
@@ -167,6 +168,15 @@ def Install(self):
         g.write('LIBS += '+self.libraries.toString(self.hwloc.lib)+'\n')
         g.write('OPTF += -DUSE_LIBHWLOC\n')
         g.write('OPTC += -DUSE_LIBHWLOC\n')
+      # To avoid a bug related to MPI_IN_PLACE and old MPICH releases, see MR 4410
+      self.avoid_mpi_in_place = 0
+      if 'download-mumps-avoid-mpi-in-place' in self.framework.clArgDB: # user-provided value takes precedence
+        self.avoid_mpi_in_place = self.clArgDB['download-mumps-avoid-mpi-in-place']
+      elif hasattr(self.mpi, 'mpich_numversion') and int(self.mpi.mpich_numversion) < 40000101:
+        self.avoid_mpi_in_place = 1
+      if self.avoid_mpi_in_place:
+        g.write('CDEFS += -DAVOID_MPI_IN_PLACE') # only take effect since mumps-5.6.2
+        self.addDefine('HAVE_MUMPS_AVOID_MPI_IN_PLACE', 1)
     g.close()
     if self.installNeeded('Makefile.inc'):
       try:

config/BuildSystem/config/packages/hip.py

@@ -8,10 +8,10 @@ class Configure(config.package.Package):
   def __init__(self, framework):
     config.package.Package.__init__(self, framework)
 
-    self.minversion       = '4.1.0'
-    # HIP_VERSION_PATCH (e.g., 31061) is not necessarily the AMD advertised patch version, e.g., in 5.6.0
-    self.versionname      = 'HIP_VERSION_MAJOR.HIP_VERSION_MINOR.HIP_VERSION_PATCH'
-    self.versioninclude   = 'hip/hip_version.h'
+    self.minversion       = '5.0.0'
+    # Check version from rocm-core here, as HIP_VERSION_PATCH (e.g., 31061 from hip_version.h) is not necessarily the AMD advertised patch version, e.g., in 5.6.0
+    self.versionname      = 'ROCM_VERSION_MAJOR.ROCM_VERSION_MINOR.ROCM_VERSION_PATCH'
+    self.versioninclude   = 'rocm_version.h'
     self.requiresversion  = 1
     self.functionsCxx     = [1,'', 'rocblas_create']
     self.includes         = ['hip/hip_runtime.h']

config/BuildSystem/config/packages/hpddm.py

@@ -4,7 +4,7 @@ class Configure(config.package.Package):
   def __init__(self,framework):
     config.package.Package.__init__(self,framework)
     self.version                = '2.2.4'
-    self.gitcommit              = '69bebeece9139b59524c70b88b60ccd4fba8a16a' # fix-saddle oct-10-2023
+    self.gitcommit              = 'f1526d323ba1dd9e1680e35a59241f6e7f89e8a1' # main oct-16-2023
     self.download               = ['git://https://github.com/hpddm/hpddm','https://github.com/hpddm/hpddm/archive/'+self.gitcommit+'.tar.gz']
     self.minversion             = '2.2.1'
     self.versionname            = 'HPDDM_VERSION'

config/BuildSystem/config/setCompilers.py

@@ -2394,7 +2394,7 @@ def checkSharedLinker(self):
         if self.sharedLibraries: break
         self.delMakeMacro('LD_SHARED')
         del self.LD_SHARED
-        del self.sharedLinker
+        if hasattr(self,'sharedLinker'): del self.sharedLinker
     return
 
   def checkLinkerFlag(self, flag, filterAlways = 0):

doc/changes/dev.rst

@@ -51,6 +51,9 @@ Changes: Development
 
 .. rubric:: TS:
 
+- Add support for custom predictor callbacks in the second-order generalized-alpha method using ``TSAlpha2SetPredictor()``
+- Allow adaptivity to change time step size in first step of second-order generalized-alpha method.
+
 .. rubric:: TAO:
 
 .. rubric:: DM/DA:

doc/community/meetings/meeting.rst

@@ -1,4 +1,3 @@
-
 .. _meetings:
 
 *********************
@@ -29,15 +28,15 @@ simulations by scientists and engineers.
 Upcoming Meetings
 =================
 
-- The 2024 user meeting will take place May 23,24 in Cologne, Germany.
+- `The 2024 user meeting will take place May 23-24 in Cologne, Germany. <https://cds.uni-koeln.de/en/workshops/petsc-2024/home>`__
 - The 2025 user meeting will take place in North America. Suggestions for a meeting location?
 
 Previous Meetings
 =================
 
-- :any:`Chicago, USA; June 5-7, 2023 <2023_meeting>`
-- `Atlanta, USA; June 5-7, 2019 <https://petsc.gitlab.io/annual-meetings/2019/index.html>`__
+- :any:`Chicago, Illinois, USA; June 5-7, 2023 <2023_meeting>`
+- `Atlanta, Georgia, USA; June 5-7, 2019 <https://petsc.gitlab.io/annual-meetings/2019/index.html>`__
 - `London, UK; June 4-6, 2018 <https://petsc.gitlab.io/annual-meetings/2018/index.html>`__
-- `Boulder, Colorado; June 14-16, 2017 <https://petsc.gitlab.io/annual-meetings/2017/index.html>`__
+- `Boulder, Colorado, USA; June 14-16, 2017 <https://petsc.gitlab.io/annual-meetings/2017/index.html>`__
 - `Vienna, Austria; June 28-30, 2016  <https://petsc.gitlab.io/annual-meetings/2016/index.html>`__
-- `Argonne National Laboratory, USA; June 15-18, 2015 <https://petsc.gitlab.io/annual-meetings/2015/index.html>`__
+- `Argonne National Laboratory, Illinois, USA; June 15-18, 2015 <https://petsc.gitlab.io/annual-meetings/2015/index.html>`__

include/petscmat.h

@@ -944,7 +944,7 @@ PETSC_EXTERN PetscErrorCode MatCreateMPIMatConcatenateSeqMat(MPI_Comm, Mat, Pets
    Note:
    For efficiency one should use `MatSetValues()` and set several values simultaneously.
 
-.seealso: [](ch_matrices), `Mat`, `MatAssemblyBegin()`, `MatAssemblyEnd`, `InsertMode`, `MatGetValue()`, `MatSetValues()`,
+.seealso: [](ch_matrices), `Mat`, `MatAssemblyBegin()`, `MatAssemblyEnd()`, `InsertMode`, `MatGetValue()`, `MatSetValues()`,
           `MatSetValueLocal()`, `MatSetValuesLocal()`
 @*/
 static inline PetscErrorCode MatSetValue(Mat mat, PetscInt i, PetscInt j, PetscScalar va, InsertMode mode)
@@ -968,13 +968,13 @@ static inline PetscErrorCode MatSetValue(Mat mat, PetscInt i, PetscInt j, PetscS
    Level: advanced
 
    Notes:
-   The matrix must have been assembled with `MatAssemblyBegin()` and `MatAssemblyEnd` before this call
+   The matrix must have been assembled with `MatAssemblyBegin()` and `MatAssemblyEnd()` before this call
 
    For efficiency one should use `MatGetValues()` and get several values simultaneously.
 
    See notes for `MatGetValues()`.
 
-.seealso: [](ch_matrices), `Mat`, `MatAssemblyBegin()`, `MatAssemblyEnd`, `MatSetValue()`, `MatGetValueLocal()`, `MatGetValues()`
+.seealso: [](ch_matrices), `Mat`, `MatAssemblyBegin()`, `MatAssemblyEnd()`, `MatSetValue()`, `MatGetValueLocal()`, `MatGetValues()`
 @*/
 static inline PetscErrorCode MatGetValue(Mat mat, PetscInt row, PetscInt col, PetscScalar *va)
 {

include/petscts.h

@@ -1371,6 +1371,25 @@ PETSC_EXTERN PetscErrorCode TSAlpha2SetRadius(TS, PetscReal);
 PETSC_EXTERN PetscErrorCode TSAlpha2SetParams(TS, PetscReal, PetscReal, PetscReal, PetscReal);
 PETSC_EXTERN PetscErrorCode TSAlpha2GetParams(TS, PetscReal *, PetscReal *, PetscReal *, PetscReal *);
 
+/*S
+  TSAlpha2Predictor - A callback to set the predictor (i.e., the initial guess for the nonlinear solver) in
+  a second-order generalized-alpha time integrator.
+
+  Calling Sequence:
+  + ts   - the `TS` context obtained from `TSCreate()`
+  . X0   - the previous time step's state vector
+  . V0   - the previous time step's first derivative of the state vector
+  . A0   - the previous time step's second derivative of the state vector
+  . X1   - the vector into which the initial guess for the current time step will be written
+  - ctx  - [optional] user-defined context for the predictor evaluation routine (may be `NULL`)
+
+  Level: intermediate
+
+  .seealso: [](ch_ts), `TS`, `TSAlpha2SetPredictor()`
+S*/
+PETSC_EXTERN_TYPEDEF typedef PetscErrorCode (*TSAlpha2Predictor)(TS ts, Vec X0, Vec V0, Vec A0, Vec X1, void *ctx);
+PETSC_EXTERN PetscErrorCode TSAlpha2SetPredictor(TS, TSAlpha2Predictor, void *ctx);
+
 PETSC_EXTERN PetscErrorCode TSSetDM(TS, DM);
 PETSC_EXTERN PetscErrorCode TSGetDM(TS, DM *);
 

src/binding/petsc4py/src/petsc4py/PETSc/Mat.pyx

@@ -90,6 +90,7 @@ class MatType(object):
     LMVMSYMBADBROYDEN = S_(MATLMVMSYMBADBROYDEN)
     LMVMDIAGBBROYDEN = S_(MATLMVMDIAGBROYDEN)
     CONSTANTDIAGONAL = S_(MATCONSTANTDIAGONAL)
+    DIAGONAL         = S_(MATDIAGONAL)
     H2OPUS           = S_(MATH2OPUS)
 
 class MatOption(object):
@@ -409,6 +410,9 @@ cdef class Mat(Object):
     def __rtruediv__(self, other):
         return mat_rdiv(self, other)
 
+    def __matmul__(self, other):
+        return mat_matmul(self, other)
+
     #
 
     def __getitem__(self, ij):
@@ -633,6 +637,25 @@ cdef class Mat(Object):
         CHKERR( MatGetVecType(self.mat, &cval) )
         return bytes2str(cval)
 
+    def setNestVecType(self, vec_type: Vec.Type | str) -> None:
+        """Set the vector type for a `Type.NEST` matrix.
+
+        Collective.
+
+        Parameters
+        ----------
+        vec_type
+            Vector type used when creating vectors with `createVecs`.
+
+        See Also
+        --------
+        getVecType, petsc.MatNestSetVecType
+
+        """
+        cdef PetscVecType cval = NULL
+        vec_type = str2bytes(vec_type, &cval)
+        CHKERR( MatNestSetVecType(self.mat, cval) )
+
     #
 
     def createAIJ(
@@ -1288,7 +1311,7 @@ cdef class Mat(Object):
 
     def createNest(
         self,
-        mats: Sequence[Mat],
+        mats: Sequence[Sequence[Mat]],
         isrows: Sequence[IS] | None = None,
         iscols: Sequence[IS] | None = None,
         comm: Comm | None = None,
@@ -1300,9 +1323,9 @@ cdef class Mat(Object):
         Parameters
         ----------
         mats
-            Row-aligned iterable of matrices with size
-            ``len(isrows)*len(iscols)``. Empty submatrices can be set with
-            `None`.
+            Iterable of matrix block rows with size ``len(isrows)``.
+            Each matrix block row must be of size ``len(iscols)``.
+            Empty submatrices can be set with `None`.
         isrows
             Index set for each nested row block, defaults to contiguous
             ordering.
@@ -1485,6 +1508,64 @@ cdef class Mat(Object):
         CHKERR( PetscCLEAR(self.obj) ); self.mat = newmat
         return self
 
+    def createConstantDiagonal(
+        self,
+        size: MatSizeSpec,
+        diag: float,
+        comm: Comm | None = None,
+        ) -> Self:
+        """Create a diagonal matrix of type `Type.CONSTANTDIAGONAL`.
+
+        Collective.
+
+        Parameters
+        ----------
+        size
+            Matrix size.
+        diag
+            The diagonal value.
+        comm
+            MPI communicator, defaults to `Sys.getDefaultComm`.
+
+        See Also
+        --------
+        createDiagonal
+
+        """
+        cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT)
+        cdef PetscInt rbs = 0, cbs = 0, m = 0, n = 0, M = 0, N = 0
+        Mat_Sizes(size, None, &rbs, &cbs, &m, &n, &M, &N)
+        Sys_Layout(ccomm, rbs, &m, &M)
+        Sys_Layout(ccomm, cbs, &n, &N)
+        cdef PetscMat newmat = NULL
+        CHKERR( MatCreateConstantDiagonal(ccomm, m, n, M, N, diag, &newmat) )
+        CHKERR( PetscCLEAR(self.obj) ); self.mat = newmat
+        return self
+
+    def createDiagonal(
+        self,
+        Vec diag,
+        ) -> Self:
+        """Create a diagonal matrix of type `Type.DIAGONAL`.
+
+        Collective.
+
+        Parameters
+        ----------
+        diag
+            The vector holding diagonal values.
+
+        See Also
+        --------
+        createConstantDiagonal
+
+        """
+        cdef PetscVec dvec = diag.vec
+        cdef PetscMat newmat = NULL
+        CHKERR( MatCreateDiagonal(dvec, &newmat) )
+        CHKERR( PetscCLEAR(self.obj) ); self.mat = newmat
+        return self
+
     def createPython(self, size: MatSizeSpec, context: Any = None, comm: Comm | None = None) -> Self:
         """Create a `Type.PYTHON` matrix.
 

src/binding/petsc4py/src/petsc4py/PETSc/SNES.pyx

@@ -1505,6 +1505,20 @@ cdef class SNES(Object):
         """
         CHKERR( SNESSetUp(self.snes) )
 
+    def setUpMatrices(self) -> None:
+        """Ensures that matrices are available for Newton-like methods.
+
+        Collective.
+
+        This is only of use to implementers of custom SNES types.
+
+        See Also
+        --------
+        setUp, petsc.SNESSetUpMatrices
+
+        """
+        CHKERR( SNESSetUpMatrices(self.snes) )
+
     def reset(self) -> None:
         """Reset the solver.
 

src/binding/petsc4py/src/petsc4py/PETSc/Vec.pyx

@@ -104,6 +104,8 @@ cdef class Vec(Object):
     def __rtruediv__(self, other):
         return vec_rdiv(self, other)
 
+    def __matmul__(self, other):
+        return vec_matmul(self, other)
     #
 
     #def __len__(self):

src/binding/petsc4py/src/petsc4py/PETSc/petscmat.pxi

@@ -83,6 +83,7 @@ cdef extern from * nogil:
     PetscMatType MATLMVMSYMBADBROYDEN
     PetscMatType MATLMVMDIAGBROYDEN
     PetscMatType MATCONSTANTDIAGONAL
+    PetscMatType MATDIAGONAL
     PetscMatType MATH2OPUS
 
     ctypedef const char* PetscMatOrderingType "MatOrderingType"
@@ -224,6 +225,8 @@ cdef extern from * nogil:
     PetscErrorCode MatCreateSeqAIJWithArrays(MPI_Comm,PetscInt,PetscInt,PetscInt[],PetscInt[],PetscScalar[],PetscMat*)
     PetscErrorCode MatCreateMPIAIJWithArrays(MPI_Comm,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt[],PetscInt[],PetscScalar[],PetscMat*)
     PetscErrorCode MatCreateMPIAIJWithSplitArrays(MPI_Comm,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt[],PetscInt[],PetscScalar[],PetscInt[],PetscInt[],PetscScalar[],PetscMat*)
+    PetscErrorCode MatCreateDiagonal(PetscVec,PetscMat*)
+    PetscErrorCode MatCreateConstantDiagonal(MPI_Comm,PetscInt,PetscInt,PetscInt,PetscInt,PetscScalar,PetscMat*)
 
     PetscErrorCode MatSetSizes(PetscMat,PetscInt,PetscInt,PetscInt,PetscInt)
     PetscErrorCode MatSetBlockSize(PetscMat,PetscInt)
@@ -271,6 +274,7 @@ cdef extern from * nogil:
     PetscErrorCode MatNestGetLocalISs(PetscMat,PetscIS*,PetscIS*)
     PetscErrorCode MatNestGetSize(PetscMat,PetscInt*,PetscInt*)
     PetscErrorCode MatNestGetSubMat(PetscMat,PetscInt,PetscInt,PetscMat*)
+    PetscErrorCode MatNestSetVecType(PetscMat,PetscVecType)
 
     PetscErrorCode MatEqual(PetscMat,PetscMat,PetscBool*)
     PetscErrorCode MatLoad(PetscMat,PetscViewer)
@@ -637,7 +641,7 @@ cdef Vec mat_mul_vec(Mat self, Vec other):
 cdef Mat mat_mul_mat(Mat self, Mat other):
     return self.matMult(other)
 
-cdef Mat mat_mul(Mat self, other):
+cdef object mat_mul(Mat self, other):
     if isinstance(other, Vec):
         return mat_mul_vec(self, <Vec>other)
     elif isinstance(other, Mat):
@@ -648,6 +652,13 @@ cdef Mat mat_mul(Mat self, other):
 cdef Mat mat_div(Mat self, other):
     return mat_idiv(mat_pos(self), other)
 
+cdef object mat_matmul(Mat self, other):
+    if isinstance(other, Vec):
+        return mat_mul_vec(self, <Vec>other)
+    if isinstance(other, Mat):
+        return mat_mul_mat(self, <Mat>other)
+    return NotImplemented
+
 # reflected binary operations
 
 cdef Mat mat_radd(Mat self, other):

src/binding/petsc4py/src/petsc4py/PETSc/petscsnes.pxi

@@ -175,6 +175,7 @@ cdef extern from * nogil:
     PetscErrorCode SNESMonitor(PetscSNES,PetscInt,PetscReal)
 
     PetscErrorCode SNESSetUp(PetscSNES)
+    PetscErrorCode SNESSetUpMatrices(PetscSNES)
     PetscErrorCode SNESReset(PetscSNES)
     PetscErrorCode SNESSolve(PetscSNES,PetscVec,PetscVec)
 

src/binding/petsc4py/src/petsc4py/PETSc/petscvec.pxi

@@ -311,6 +311,15 @@ cdef Vec vec_mul(Vec self, other):
 cdef Vec vec_div(Vec self, other):
     return vec_idiv(vec_pos(self), other)
 
+cdef object vec_matmul(Vec self, other):
+    if isinstance(other, Vec):
+        return self.dot(other)
+    if isinstance(other, Mat):
+        result = other.createVecRight()
+        other.multTranspose(self, result)
+        return result
+    return NotImplemented
+
 # reflected binary operations
 
 cdef Vec vec_radd(Vec self, other):