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