MUMPS FetchContent CMake use example

Readme.md

MUMPS example via CMake FetchContent

Example of using MUMPS via CMake FetchContent.

This example will download and build MUMPS, and find the required Scalapack and LAPACK. If Scalapack and Lapack aren't available, they can be built via our MUMPS CMake project scripts/CMakeLists.txt.

It requires that MPI library for Fortran and C is already installed.

cmake -B build

cmake --build build

# optional
ctest --test-dir build

CMakeLists.txt

cmake_minimum_required(VERSION 3.14...3.26)

project(MUMPS_example
LANGUAGES C Fortran
)

enable_testing()

include(FetchContent)

set(FETCHCONTENT_UPDATES_DISCONNECTED_MUMPSUPSTREAM true)
set(FETCHCONTENT_QUIET false)

file(GENERATE OUTPUT ${PROJECT_BINARY_DIR}/.gitignore CONTENT "*")

# --- build MUMPS

find_package(MPI COMPONENTS C Fortran REQUIRED)

FetchContent_Declare(MUMPSupstream
  GIT_REPOSITORY https://github.com/scivision/mumps.git
  GIT_TAG v5.5.1.11
  TLS_VERIFY on
  INACTIVITY_TIMEOUT 60
)

FetchContent_MakeAvailable(MUMPSupstream)

# Example C and Fortran main programs using MUMPS

add_executable(c_mumps d_simple.c)
target_link_libraries(c_mumps PRIVATE MUMPS::MUMPS MPI::MPI_C)
add_test(NAME MUMPS_C
COMMAND ${MPIEXEC_EXECUTABLE} ${MPIEXEC_NUMPROC_FLAG} 2 $<TARGET_FILE:c_mumps>
)

add_executable(fortran_mumps d_simple.f90)
target_link_libraries(fortran_mumps PRIVATE MUMPS::MUMPS MPI::MPI_Fortran)
add_test(NAME MUMPS_Fortran
COMMAND ${MPIEXEC_EXECUTABLE} ${MPIEXEC_NUMPROC_FLAG} 2 $<TARGET_FILE:fortran_mumps>
)

d_simple.c

/*
 *
 *  This file is part of MUMPS 5.2.1, released
 *  on Fri Jun 14 14:46:05 UTC 2019
 *
 */
/* Example program using the C interface to the
 * double real arithmetic version of MUMPS, dmumps_c.
 * We solve the system A x = RHS with
 *   A = diag(1 2) and RHS = [1 4]^T
 * Solution is [1 2]^T */
#include <stdio.h>
#include <string.h>

#include "mpi.h"
#include "dmumps_c.h"

#define JOB_INIT -1
#define JOB_END -2
#define USE_COMM_WORLD -987654


int main(int argc, char ** argv)
{
  DMUMPS_STRUC_C id;
  MUMPS_INT n = 2;
  MUMPS_INT8 nnz = 2;
  MUMPS_INT irn[] = {1,2};
  MUMPS_INT jcn[] = {1,2};
  double a[2];
  double rhs[2];

  MUMPS_INT myid, ierr;

  ierr = MPI_Init(&argc, &argv);
  if (ierr != 0) {
     fprintf(stderr, "failed to init MPI");
     return 1;
  }
  ierr = MPI_Comm_rank(MPI_COMM_WORLD, &myid);
  /* Define A and rhs */
  rhs[0]=1.0;rhs[1]=4.0;
  a[0]=1.0;a[1]=2.0;

  /* Initialize a MUMPS instance. Use MPI_COMM_WORLD */
  id.comm_fortran=USE_COMM_WORLD;
  id.par=1; id.sym=0;
  id.job=JOB_INIT;
  dmumps_c(&id);

  /* Define the problem on the host */
  if (myid == 0) {
    id.n = n; id.nnz =nnz; id.irn=irn; id.jcn=jcn;
    id.a = a; id.rhs = rhs;
  }
#define ICNTL(I) icntl[(I)-1] /* macro s.t. indices match documentation */
  /* No outputs */
  id.ICNTL(1)=-1; id.ICNTL(2)=-1; id.ICNTL(3)=-1; id.ICNTL(4)=0;

  /* Call the MUMPS package (analyse, factorization and solve). */
  id.job=6;
  dmumps_c(&id);
  if (id.infog[0]<0) {
    fprintf(stderr," (PROC %d) ERROR RETURN: \tINFOG(1)= %d\n\t\t\t\tINFOG(2)= %d\n",
        myid, id.infog[0], id.infog[1]);
    return 1;
  }

  /* Terminate instance. */
  id.job=JOB_END;
  dmumps_c(&id);
  if (myid == 0) printf("Solution is : (%8.2f  %8.2f)\n", rhs[0],rhs[1]);
  ierr = MPI_Finalize();
  return 0;
}

d_simple.f90

program d_simple
!  This file is part of MUMPS 5.2.1, released
!  on Fri Jun 14 14:46:05 UTC 2019

use, intrinsic :: iso_fortran_env, only: stderr=>error_unit, stdout=>output_unit, int64

use mpi

implicit none

INCLUDE 'dmumps_struc.h'
TYPE (DMUMPS_STRUC) mumps_par
INTEGER :: IERR, I, num_mpi
INTEGER(int64) :: I8

CALL MPI_INIT(IERR)
if(ierr/=0) error stop 'problem initializing MPI'

call MPI_COMM_size(MPI_COMM_WORLD, num_mpi, ierr)
if(ierr/=0) error stop 'problem getting number of MPI processes'
print '(A,I3,A)', 'using ',num_mpi,' MPI processes'
! Define a communicator for the package.
mumps_par%COMM = MPI_COMM_WORLD
!  Initialize an instance of the package
!  for L U factorization (sym = 0, with working host)
mumps_par%JOB = -1
mumps_par%SYM = 0
mumps_par%PAR = 1

CALL DMUMPS(mumps_par)

mumps_par%icntl(1) = stderr  ! error messages
mumps_par%icntl(2) = stdout !  diagnosic, statistics, and warning messages
mumps_par%icntl(3) = stdout! ! global info, for the host (myid==0)
mumps_par%icntl(4) = 1           ! default is 2, this reduces verbosity

! === config done, now check config
IF (mumps_par%INFOG(1) < 0) THEN
  WRITE(stderr,'(A,A,I6,A,I9)') " ERROR RETURN: ", &
  "  mumps_par%INFOG(1)= ", mumps_par%INFOG(1), &
  "  mumps_par%INFOG(2)= ", mumps_par%INFOG(2)

  error stop
END IF
!  Define problem on the host (processor 0)
IF ( mumps_par%MYID == 0 ) THEN
  mumps_par%N = 5
  mumps_par%NNZ = 12
  ALLOCATE( mumps_par%IRN ( mumps_par%NNZ ) )
  ALLOCATE( mumps_par%JCN ( mumps_par%NNZ ) )
  ALLOCATE( mumps_par%A( mumps_par%NNZ ) )
  ALLOCATE( mumps_par%RHS ( mumps_par%N  ) )
  mumps_par%IRN = [1,2,4,5,2,1,5,3,2,3,1,3]
  mumps_par%JCN = [2,3,3,5,1,1,2,4,5,2,3,3]
  mumps_par%A   = [3., -3., 2., 1., 3., 2., 4., 2., 6., -1., 4., 1.]

  mumps_par%RHS = [20., 24., 9., 6., 13.]
END IF
!  Call package for solution
mumps_par%JOB = 6
CALL DMUMPS(mumps_par)
IF (mumps_par%INFOG(1) < 0) THEN
  WRITE(stderr,'(A,A,I6,A,I9)') " ERROR RETURN: ", &
  "  mumps_par%INFOG(1)= ", mumps_par%INFOG(1), &
  "  mumps_par%INFOG(2)= ", mumps_par%INFOG(2)

 error stop
END IF
!  Solution has been assembled on the host

IF ( mumps_par%MYID == 0 ) THEN
  print *, ' Solution is: '
  print *, mumps_par%RHS

  if (sum(mumps_par%rhs-[1,2,3,4,5]) > 0.01) error stop 'excessive error in computation'
END IF
!  Deallocate user data
IF ( mumps_par%MYID == 0 )THEN
  DEALLOCATE( mumps_par%IRN )
  DEALLOCATE( mumps_par%JCN )
  DEALLOCATE( mumps_par%A   )
  DEALLOCATE( mumps_par%RHS )
END IF
!  Destroy the instance (deallocate internal data structures)
mumps_par%JOB = -2
CALL DMUMPS(mumps_par)
IF (mumps_par%INFOG(1) < 0) THEN
  WRITE(stderr,'(A,A,I6,A,I9)') " ERROR RETURN: ", &
  "  mumps_par%INFOG(1)= ", mumps_par%INFOG(1), &
  "  mumps_par%INFOG(2)= ", mumps_par%INFOG(2)

 error stop
END IF

call mpi_finalize(ierr)

END program