1 2 /* 3 Code for some of the parallel vector primatives. 4 */ 5 #include <../src/vec/vec/impls/mpi/pvecimpl.h> 6 #include <petscblaslapack.h> 7 8 PetscErrorCode VecMDot_MPI(Vec xin,PetscInt nv,const Vec y[],PetscScalar *z) 9 { 10 PetscScalar awork[128],*work = awork; 11 PetscErrorCode ierr; 12 13 PetscFunctionBegin; 14 if (nv > 128) { 15 ierr = PetscMalloc1(nv,&work);CHKERRQ(ierr); 16 } 17 ierr = VecMDot_Seq(xin,nv,y,work);CHKERRQ(ierr); 18 ierr = MPIU_Allreduce(work,z,nv,MPIU_SCALAR,MPIU_SUM,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 19 if (nv > 128) { 20 ierr = PetscFree(work);CHKERRQ(ierr); 21 } 22 PetscFunctionReturn(0); 23 } 24 25 PetscErrorCode VecMTDot_MPI(Vec xin,PetscInt nv,const Vec y[],PetscScalar *z) 26 { 27 PetscScalar awork[128],*work = awork; 28 PetscErrorCode ierr; 29 30 PetscFunctionBegin; 31 if (nv > 128) { 32 ierr = PetscMalloc1(nv,&work);CHKERRQ(ierr); 33 } 34 ierr = VecMTDot_Seq(xin,nv,y,work);CHKERRQ(ierr); 35 ierr = MPIU_Allreduce(work,z,nv,MPIU_SCALAR,MPIU_SUM,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 36 if (nv > 128) { 37 ierr = PetscFree(work);CHKERRQ(ierr); 38 } 39 PetscFunctionReturn(0); 40 } 41 42 #include <../src/vec/vec/impls/seq/ftn-kernels/fnorm.h> 43 PetscErrorCode VecNorm_MPI(Vec xin,NormType type,PetscReal *z) 44 { 45 PetscReal sum,work = 0.0; 46 const PetscScalar *xx; 47 PetscErrorCode ierr; 48 PetscInt n = xin->map->n; 49 PetscBLASInt one = 1,bn = 0; 50 51 PetscFunctionBegin; 52 ierr = PetscBLASIntCast(n,&bn);CHKERRQ(ierr); 53 if (type == NORM_2 || type == NORM_FROBENIUS) { 54 ierr = VecGetArrayRead(xin,&xx);CHKERRQ(ierr); 55 work = PetscRealPart(BLASdot_(&bn,xx,&one,xx,&one)); 56 ierr = VecRestoreArrayRead(xin,&xx);CHKERRQ(ierr); 57 ierr = MPIU_Allreduce(&work,&sum,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 58 *z = PetscSqrtReal(sum); 59 ierr = PetscLogFlops(2.0*xin->map->n);CHKERRQ(ierr); 60 } else if (type == NORM_1) { 61 /* Find the local part */ 62 ierr = VecNorm_Seq(xin,NORM_1,&work);CHKERRQ(ierr); 63 /* Find the global max */ 64 ierr = MPIU_Allreduce(&work,z,1,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 65 } else if (type == NORM_INFINITY) { 66 /* Find the local max */ 67 ierr = VecNorm_Seq(xin,NORM_INFINITY,&work);CHKERRQ(ierr); 68 /* Find the global max */ 69 ierr = MPIU_Allreduce(&work,z,1,MPIU_REAL,MPIU_MAX,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 70 } else if (type == NORM_1_AND_2) { 71 PetscReal temp[2]; 72 ierr = VecNorm_Seq(xin,NORM_1,temp);CHKERRQ(ierr); 73 ierr = VecNorm_Seq(xin,NORM_2,temp+1);CHKERRQ(ierr); 74 temp[1] = temp[1]*temp[1]; 75 ierr = MPIU_Allreduce(temp,z,2,MPIU_REAL,MPIU_SUM,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 76 z[1] = PetscSqrtReal(z[1]); 77 } 78 PetscFunctionReturn(0); 79 } 80 81 PetscErrorCode VecMax_MPI(Vec xin,PetscInt *idx,PetscReal *z) 82 { 83 PetscErrorCode ierr; 84 PetscReal work; 85 86 PetscFunctionBegin; 87 /* Find the local max */ 88 ierr = VecMax_Seq(xin,idx,&work);CHKERRQ(ierr); 89 90 /* Find the global max */ 91 if (!idx) { 92 ierr = MPIU_Allreduce(&work,z,1,MPIU_REAL,MPIU_MAX,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 93 } else { 94 PetscReal work2[2],z2[2]; 95 PetscInt rstart; 96 rstart = xin->map->rstart; 97 work2[0] = work; 98 work2[1] = *idx + rstart; 99 ierr = MPIU_Allreduce(work2,z2,2,MPIU_REAL,MPIU_MAXINDEX_OP,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 100 *z = z2[0]; 101 *idx = (PetscInt)z2[1]; 102 } 103 PetscFunctionReturn(0); 104 } 105 106 PetscErrorCode VecMin_MPI(Vec xin,PetscInt *idx,PetscReal *z) 107 { 108 PetscErrorCode ierr; 109 PetscReal work; 110 111 PetscFunctionBegin; 112 /* Find the local Min */ 113 ierr = VecMin_Seq(xin,idx,&work);CHKERRQ(ierr); 114 115 /* Find the global Min */ 116 if (!idx) { 117 ierr = MPIU_Allreduce(&work,z,1,MPIU_REAL,MPIU_MIN,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 118 } else { 119 PetscReal work2[2],z2[2]; 120 PetscInt rstart; 121 122 ierr = VecGetOwnershipRange(xin,&rstart,NULL);CHKERRQ(ierr); 123 work2[0] = work; 124 work2[1] = *idx + rstart; 125 ierr = MPIU_Allreduce(work2,z2,2,MPIU_REAL,MPIU_MININDEX_OP,PetscObjectComm((PetscObject)xin));CHKERRQ(ierr); 126 *z = z2[0]; 127 *idx = (PetscInt)z2[1]; 128 } 129 PetscFunctionReturn(0); 130 } 131