Actual source code: ex54f.F90
1: ! test verifies DMShellSetCreateFieldDecomposition interface in Fortran
2: program main
3: #include "petsc/finclude/petsc.h"
4: #include "petsc/finclude/petscdmshell.h"
6: use petsc
7: use petscdmshell
8: implicit none
9: type(tDM) :: dm
10: PetscErrorCode :: ierr
11: interface
12: subroutine myFieldDecomp(dm, nfields, fieldNames, isFields, subDms, ierr)
13: use petsc
14: implicit none
15: type(tDM), intent(in) :: dm
16: PetscInt, intent(out) :: nfields
17: character(len=30), allocatable, intent(out) :: fieldNames(:)
18: type(tIS), allocatable, intent(out) :: isFields(:)
19: type(tDM), allocatable, intent(out) :: subDms(:)
20: PetscErrorCode, intent(out) :: ierr
21: end subroutine myFieldDecomp
22: end interface
23: ! initializing PETSc
24: PetscCallA(PetscInitialize(PETSC_NULL_CHARACTER, ierr))
25: ! creating a DMShell object
26: PetscCallA(DMShellCreate(PETSC_COMM_WORLD, dm, ierr))
27: ! registering the Fortran field decomposition callback
28: PetscCallA(DMShellSetCreateFieldDecomposition(dm, myFieldDecomp, ierr))
29: ! for this minimal test, we simply print a success message to the console
30: print *, 'DMShellSetCreateFieldDecomposition set successfully.'
31: ! cleanup
32: PetscCallA(DMDestroy(dm, ierr))
33: PetscCallA(PetscFinalize(ierr))
34: end program main
36: ! a simple Fortran callback for field decomposition.
37: subroutine myFieldDecomp(dm, nfields, fieldNames, isFields, subDms, ierr)
38: use petsc
39: implicit none
40: type(tDM), intent(in) :: dm
41: PetscInt, intent(out) :: nfields
42: character(len=30), allocatable, intent(out) :: fieldNames(:)
43: type(tIS), allocatable, intent(out) :: isFields(:)
44: type(tDM), allocatable, intent(out) :: subDms(:)
45: PetscErrorCode, intent(out) :: ierr
46: PetscInt :: i
47: ! defining a simple decomposition with two fields
48: nfields = 2
49: allocate (fieldNames(nfields))
50: allocate (isFields(nfields))
51: allocate (subDms(nfields))
52: fieldNames(1) = 'field1'
53: fieldNames(2) = 'field2'
54: ! set the pointer arrays to NULL (using pointer assignment)
55: do i = 1, nfields
56: isFields(i) = PETSC_NULL_IS
57: subDms(i) = PETSC_NULL_DM
58: end do
59: ierr = 0
60: print *, 'myFieldDecomp callback invoked.'
61: end subroutine myFieldDecomp
62: !/*TEST
63: !
64: ! test:
65: !TEST*/