1 program main 2!----------------------------------------------------------------------- 3! 4! Tests DMDAGetVecGetArray() 5!----------------------------------------------------------------------- 6! 7 8#include <petsc/finclude/petscdm.h> 9 use petsc 10 implicit none 11 12 Type(tVec) g 13 Type(tDM) ada 14 15 PetscScalar,pointer :: x1(:),x2(:,:) 16 PetscScalar,pointer :: x3(:,:,:),x4(:,:,:,:) 17 PetscErrorCode ierr 18 PetscInt m,n,p,dof,s,i,j,k,xs,xl 19 PetscInt ys,yl 20 PetscInt zs,zl,sw 21 22 PetscInt nen,nel 23 PetscInt, pointer :: elements(:) 24 25 m = 5 26 n = 6 27 p = 4; 28 s = 1 29 dof = 1 30 sw = 1 31 PetscCallA(PetscInitialize(ierr)) 32 PetscCallA(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m,dof,sw,PETSC_NULL_INTEGER,ada,ierr)) 33 PetscCallA(DMSetUp(ada,ierr)) 34 PetscCallA(DMGetGlobalVector(ada,g,ierr)) 35 PetscCallA(DMDAGetCorners(ada,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,xl,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr)) 36 PetscCallA(DMDAVecGetArrayF90(ada,g,x1,ierr)) 37 do i=xs,xs+xl-1 38! CHKMEMQ 39 x1(i) = i 40! CHKMEMQ 41 enddo 42 PetscCallA(DMDAVecRestoreArrayF90(ada,g,x1,ierr)) 43 PetscCallA(VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)) 44 PetscCallA(DMRestoreGlobalVector(ada,g,ierr)) 45 PetscCallA(DMDestroy(ada,ierr)) 46 47 PetscCallA(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,m,n,PETSC_DECIDE,PETSC_DECIDE,dof,s,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr)) 48 PetscCallA(DMSetUp(ada,ierr)) 49 PetscCallA(DMGetGlobalVector(ada,g,ierr)) 50 PetscCallA(DMDAGetCorners(ada,xs,ys,PETSC_NULL_INTEGER,xl,yl,PETSC_NULL_INTEGER,ierr)) 51 PetscCallA(DMDAVecGetArrayF90(ada,g,x2,ierr)) 52 do i=xs,xs+xl-1 53 do j=ys,ys+yl-1 54! CHKMEMQ 55 x2(i,j) = i + j 56! CHKMEMQ 57 enddo 58 enddo 59 PetscCallA(DMDAVecRestoreArrayF90(ada,g,x2,ierr)) 60 PetscCallA(VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)) 61 PetscCallA(DMRestoreGlobalVector(ada,g,ierr)) 62 63 PetscCallA(DMDAGetElements(ada,nen,nel,elements,ierr)) 64 do i=1,nen*nel 65 PetscCheckA(elements(i) .ge. 0,PETSC_COMM_SELF,PETSC_ERR_PLIB,'Error getting DMDA elements') 66 enddo 67 PetscCallA(DMDARestoreElements(ada,nen,nel,elements,ierr)) 68 PetscCallA(DMDestroy(ada,ierr)) 69 70 PetscCallA(DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX, m,n,p,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,s,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr)) 71 PetscCallA(DMSetUp(ada,ierr)) 72 PetscCallA(DMGetGlobalVector(ada,g,ierr)) 73 PetscCallA(DMDAGetCorners(ada,xs,ys,zs,xl,yl,zl,ierr)) 74 PetscCallA(DMDAVecGetArrayF90(ada,g,x3,ierr)) 75 do i=xs,xs+xl-1 76 do j=ys,ys+yl-1 77 do k=zs,zs+zl-1 78! CHKMEMQ 79 x3(i,j,k) = i + j + k 80! CHKMEMQ 81 enddo 82 enddo 83 enddo 84 PetscCallA(DMDAVecRestoreArrayF90(ada,g,x3,ierr)) 85 PetscCallA(VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)) 86 PetscCallA(DMRestoreGlobalVector(ada,g,ierr)) 87 PetscCallA(DMDestroy(ada,ierr)) 88 89! 90! Same tests but now with DOF > 1, so dimensions of array are one higher 91! 92 dof = 2 93 PetscCallA(DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m,dof,sw,PETSC_NULL_INTEGER,ada,ierr)) 94 PetscCallA(DMSetUp(ada,ierr)) 95 PetscCallA(DMGetGlobalVector(ada,g,ierr)) 96 PetscCallA(DMDAGetCorners(ada,xs,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,xl,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ierr)) 97 PetscCallA(DMDAVecGetArrayF90(ada,g,x2,ierr)) 98 do i=xs,xs+xl-1 99! CHKMEMQ 100 x2(0,i) = i 101 x2(1,i) = -i 102! CHKMEMQ 103 enddo 104 PetscCallA(DMDAVecRestoreArrayF90(ada,g,x1,ierr)) 105 PetscCallA(VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)) 106 PetscCallA(DMRestoreGlobalVector(ada,g,ierr)) 107 PetscCallA(DMDestroy(ada,ierr)) 108 109 dof = 2 110 PetscCallA(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,m,n,PETSC_DECIDE,PETSC_DECIDE,dof,s,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr)) 111 PetscCallA(DMSetUp(ada,ierr)) 112 PetscCallA(DMGetGlobalVector(ada,g,ierr)) 113 PetscCallA(DMDAGetCorners(ada,xs,ys,PETSC_NULL_INTEGER,xl,yl,PETSC_NULL_INTEGER,ierr)) 114 PetscCallA(DMDAVecGetArrayF90(ada,g,x3,ierr)) 115 do i=xs,xs+xl-1 116 do j=ys,ys+yl-1 117! CHKMEMQ 118 x3(0,i,j) = i + j 119 x3(1,i,j) = -(i + j) 120! CHKMEMQ 121 enddo 122 enddo 123 PetscCallA(DMDAVecRestoreArrayF90(ada,g,x3,ierr)) 124 PetscCallA(VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)) 125 PetscCallA(DMRestoreGlobalVector(ada,g,ierr)) 126 PetscCallA(DMDestroy(ada,ierr)) 127 128 dof = 3 129 PetscCallA(DMDACreate3d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_BOX,m,n,p,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,s,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,PETSC_NULL_INTEGER,ada,ierr)) 130 PetscCallA(DMSetUp(ada,ierr)) 131 PetscCallA(DMGetGlobalVector(ada,g,ierr)) 132 PetscCallA(DMDAGetCorners(ada,xs,ys,zs,xl,yl,zl,ierr)) 133 PetscCallA(DMDAVecGetArrayF90(ada,g,x4,ierr)) 134 do i=xs,xs+xl-1 135 do j=ys,ys+yl-1 136 do k=zs,zs+zl-1 137! CHKMEMQ 138 x4(0,i,j,k) = i + j + k 139 x4(1,i,j,k) = -(i + j + k) 140 x4(2,i,j,k) = i + j + k 141! CHKMEMQ 142 enddo 143 enddo 144 enddo 145 PetscCallA(DMDAVecRestoreArrayF90(ada,g,x4,ierr)) 146 PetscCallA(VecView(g,PETSC_VIEWER_STDOUT_WORLD,ierr)) 147 PetscCallA(DMRestoreGlobalVector(ada,g,ierr)) 148 PetscCallA(DMDestroy(ada,ierr)) 149 150 PetscCallA(PetscFinalize(ierr)) 151 END PROGRAM 152 153! 154!/*TEST 155! 156! build: 157! requires: !complex 158! 159! test: 160! filter: Error: grep -v "Vec Object" | grep -v "Warning: ieee_inexact is signaling" 161! 162!TEST*/ 163