Actual source code: ex49f.F90

  1: !
  2: !  Test Fortran binding of sort routines
  3: !
  4: module ex49fmodule
  5:   use petsc
  6: #include "petsc/finclude/petsc.h"
  7:   implicit none
  8:   type uctx
  9:      PetscInt myint
 10:   end type uctx
 11: contains
 12:   subroutine CompareIntegers(a,b,ctx,res)
 13:     implicit none

 15:     PetscInt :: a,b
 16:     type(uctx) :: ctx
 17:     integer  :: res

 19:     if (a < b) then
 20:        res = -1
 21:     else if (a == b) then
 22:        res = 0
 23:     else
 24:        res = 1
 25:     end if
 26:   end subroutine CompareIntegers
 27: end module ex49fmodule

 29: program main

 31:   use ex49fmodule
 32:   implicit none

 34:   PetscErrorCode          ierr
 35:   PetscCount,parameter::  iN = 3
 36:   PetscInt, parameter ::  N = 3
 37:   PetscInt                x(N),x1(N),y(N),z(N)
 38:   PetscMPIInt             mx(N),my(N)
 39:   PetscScalar             s(N)
 40:   PetscReal               r(N)
 41:   PetscMPIInt,parameter:: two=2, five=5, seven=7
 42:   type(uctx)::            ctx
 43:   PetscInt                i
 44:   PetscSizeT              sizeofentry

 46:   PetscCallA(PetscInitialize(ierr))

 48:   x  = [3, 2, 1]
 49:   x1 = [3, 2, 1]
 50:   y  = [6, 5, 4]
 51:   z  = [3, 5, 2]
 52:   mx = [five, seven, two]
 53:   my = [five, seven, two]
 54:   s  = [1.0, 2.0, 3.0]
 55:   r  = [1.0, 2.0, 3.0]
 56: #if defined(PETSC_USE_64BIT_INDICES)
 57:   sizeofentry = 8;
 58: #else
 59:   sizeofentry = 4;
 60: #endif
 61:   ctx%myint = 1
 62:   PetscCallA(PetscSortInt(iN,x,ierr))
 63:   PetscCallA(PetscTimSort(N,x1,sizeofentry,CompareIntegers,ctx,ierr))
 64:   do i = 1,N
 65:     PetscCheckA(x1(i) .eq. x(i),PETSC_COMM_SELF,PETSC_ERR_PLIB,'PetscTimSort and PetscSortInt arrays did not match')
 66:   end do
 67:   PetscCallA(PetscSortIntWithArray(iN,y,x,ierr))
 68:   PetscCallA(PetscSortIntWithArrayPair(iN,x,y,z,ierr))

 70:   PetscCallA(PetscSortMPIInt(iN,mx,ierr))
 71:   PetscCallA(PetscSortMPIIntWithArray(iN,mx,my,ierr))
 72:   PetscCallA(PetscSortMPIIntWithIntArray(iN,mx,y,ierr))

 74:   PetscCallA(PetscSortIntWithScalarArray(iN,x,s,ierr))

 76:   PetscCallA(PetscSortReal(iN,r,ierr))
 77:   PetscCallA(PetscSortRealWithArrayInt(iN,r,x,ierr))

 79:   PetscCallA(PetscFinalize(ierr))
 80: end program main

 82: !/*TEST
 83: !
 84: !   test:
 85: !
 86: !TEST*/