1 #include <petsc/private/fortranimpl.h> 2 3 /*MC 4 PetscFortranAddr - a variable type in Fortran that can hold a 5 regular C pointer. 6 7 Notes: 8 Used, for example, as the file argument in PetscFOpen() 9 10 Level: beginner 11 12 .seealso: PetscOffset, PetscInt 13 M*/ 14 /*MC 15 PetscOffset - a variable type in Fortran used with VecGetArray() 16 and ISGetIndices() 17 18 Level: beginner 19 20 .seealso: PetscFortranAddr, PetscInt 21 M*/ 22 23 /* 24 This is code for translating PETSc memory addresses to integer offsets 25 for Fortran. 26 */ 27 char *PETSC_NULL_CHARACTER_Fortran = 0; 28 void *PETSC_NULL_INTEGER_Fortran = 0; 29 void *PETSC_NULL_SCALAR_Fortran = 0; 30 void *PETSC_NULL_DOUBLE_Fortran = 0; 31 void *PETSC_NULL_REAL_Fortran = 0; 32 void *PETSC_NULL_BOOL_Fortran = 0; 33 EXTERN_C_BEGIN 34 void (*PETSC_NULL_FUNCTION_Fortran)(void) = 0; 35 EXTERN_C_END 36 37 size_t PetscIntAddressToFortran(const PetscInt *base,const PetscInt *addr) 38 { 39 size_t tmp1 = (size_t) base,tmp2 = 0; 40 size_t tmp3 = (size_t) addr; 41 size_t itmp2; 42 43 #if !defined(PETSC_HAVE_CRAY90_POINTER) 44 if (tmp3 > tmp1) { 45 tmp2 = (tmp3 - tmp1)/sizeof(PetscInt); 46 itmp2 = (size_t) tmp2; 47 } else { 48 tmp2 = (tmp1 - tmp3)/sizeof(PetscInt); 49 itmp2 = -((size_t) tmp2); 50 } 51 #else 52 if (tmp3 > tmp1) { 53 tmp2 = (tmp3 - tmp1); 54 itmp2 = (size_t) tmp2; 55 } else { 56 tmp2 = (tmp1 - tmp3); 57 itmp2 = -((size_t) tmp2); 58 } 59 #endif 60 61 if (base + itmp2 != addr) { 62 (*PetscErrorPrintf)("PetscIntAddressToFortran:C and Fortran arrays are\n"); 63 (*PetscErrorPrintf)("not commonly aligned or are too far apart to be indexed \n"); 64 (*PetscErrorPrintf)("by an integer. Locations: C %uld Fortran %uld\n",tmp1,tmp3); 65 PETSCABORT(PETSC_COMM_WORLD,PETSC_ERR_PLIB); 66 } 67 return itmp2; 68 } 69 70 PetscInt *PetscIntAddressFromFortran(const PetscInt *base,size_t addr) 71 { 72 return (PetscInt *)(base + addr); 73 } 74 75 /* 76 obj - PETSc object on which request is made 77 base - Fortran array address 78 addr - C array address 79 res - will contain offset from C to Fortran 80 shift - number of bytes that prevent base and addr from being commonly aligned 81 N - size of the array 82 83 align indicates alignment relative to PetscScalar, 1 means aligned on PetscScalar, 2 means aligned on 2 PetscScalar etc 84 */ 85 PetscErrorCode PetscScalarAddressToFortran(PetscObject obj,PetscInt align,PetscScalar *base,PetscScalar *addr,PetscInt N,size_t *res) 86 { 87 size_t tmp1 = (size_t) base,tmp2; 88 size_t tmp3 = (size_t) addr; 89 size_t itmp2; 90 PetscInt shift; 91 92 #if !defined(PETSC_HAVE_CRAY90_POINTER) 93 if (tmp3 > tmp1) { /* C is bigger than Fortran */ 94 tmp2 = (tmp3 - tmp1)/sizeof(PetscScalar); 95 itmp2 = (size_t) tmp2; 96 shift = (align*sizeof(PetscScalar) - (PetscInt)((tmp3 - tmp1) % (align*sizeof(PetscScalar)))) % (align*sizeof(PetscScalar)); 97 } else { 98 tmp2 = (tmp1 - tmp3)/sizeof(PetscScalar); 99 itmp2 = -((size_t) tmp2); 100 shift = (PetscInt)((tmp1 - tmp3) % (align*sizeof(PetscScalar))); 101 } 102 #else 103 if (tmp3 > tmp1) { /* C is bigger than Fortran */ 104 tmp2 = (tmp3 - tmp1); 105 itmp2 = (size_t) tmp2; 106 } else { 107 tmp2 = (tmp1 - tmp3); 108 itmp2 = -((size_t) tmp2); 109 } 110 shift = 0; 111 #endif 112 113 if (shift) { 114 /* 115 Fortran and C not PetscScalar aligned,recover by copying values into 116 memory that is aligned with the Fortran 117 */ 118 PetscErrorCode ierr; 119 PetscScalar *work; 120 PetscContainer container; 121 122 ierr = PetscMalloc1(N+align,&work);CHKERRQ(ierr); 123 124 /* recompute shift for newly allocated space */ 125 tmp3 = (size_t) work; 126 if (tmp3 > tmp1) { /* C is bigger than Fortran */ 127 shift = (align*sizeof(PetscScalar) - (PetscInt)((tmp3 - tmp1) % (align*sizeof(PetscScalar)))) % (align*sizeof(PetscScalar)); 128 } else { 129 shift = (PetscInt)((tmp1 - tmp3) % (align*sizeof(PetscScalar))); 130 } 131 132 /* shift work by that number of bytes */ 133 work = (PetscScalar*)(((char*)work) + shift); 134 ierr = PetscArraycpy(work,addr,N);CHKERRQ(ierr); 135 136 /* store in the first location in addr how much you shift it */ 137 ((PetscInt*)addr)[0] = shift; 138 139 ierr = PetscContainerCreate(PETSC_COMM_SELF,&container);CHKERRQ(ierr); 140 ierr = PetscContainerSetPointer(container,addr);CHKERRQ(ierr); 141 ierr = PetscObjectCompose(obj,"GetArrayPtr",(PetscObject)container);CHKERRQ(ierr); 142 143 tmp3 = (size_t) work; 144 if (tmp3 > tmp1) { /* C is bigger than Fortran */ 145 tmp2 = (tmp3 - tmp1)/sizeof(PetscScalar); 146 itmp2 = (size_t) tmp2; 147 shift = (align*sizeof(PetscScalar) - (PetscInt)((tmp3 - tmp1) % (align*sizeof(PetscScalar)))) % (align*sizeof(PetscScalar)); 148 } else { 149 tmp2 = (tmp1 - tmp3)/sizeof(PetscScalar); 150 itmp2 = -((size_t) tmp2); 151 shift = (PetscInt)((tmp1 - tmp3) % (align*sizeof(PetscScalar))); 152 } 153 if (shift) { 154 (*PetscErrorPrintf)("PetscScalarAddressToFortran:C and Fortran arrays are\n"); 155 (*PetscErrorPrintf)("not commonly aligned.\n"); 156 (*PetscErrorPrintf)("Locations/sizeof(PetscScalar): C %f Fortran %f\n",((PetscReal)tmp3)/(PetscReal)sizeof(PetscScalar),((PetscReal)tmp1)/(PetscReal)sizeof(PetscScalar)); 157 PETSCABORT(PETSC_COMM_WORLD,PETSC_ERR_PLIB); 158 } 159 ierr = PetscInfo(obj,"Efficiency warning, copying array in XXXGetArray() due\n\ 160 to alignment differences between C and Fortran\n");CHKERRQ(ierr); 161 } 162 *res = itmp2; 163 return 0; 164 } 165 166 /* 167 obj - the PETSc object where the scalar pointer came from 168 base - the Fortran array address 169 addr - the Fortran offset from base 170 N - the amount of data 171 172 lx - the array space that is to be passed to XXXXRestoreArray() 173 */ 174 PetscErrorCode PetscScalarAddressFromFortran(PetscObject obj,PetscScalar *base,size_t addr,PetscInt N,PetscScalar **lx) 175 { 176 PetscErrorCode ierr; 177 PetscInt shift; 178 PetscContainer container; 179 PetscScalar *tlx; 180 181 ierr = PetscObjectQuery(obj,"GetArrayPtr",(PetscObject*)&container);CHKERRQ(ierr); 182 if (container) { 183 ierr = PetscContainerGetPointer(container,(void**)lx);CHKERRQ(ierr); 184 tlx = base + addr; 185 186 shift = *(PetscInt*)*lx; 187 ierr = PetscArraycpy(*lx,tlx,N);CHKERRQ(ierr); 188 tlx = (PetscScalar*)(((char*)tlx) - shift); 189 190 ierr = PetscFree(tlx);CHKERRQ(ierr); 191 ierr = PetscContainerDestroy(&container);CHKERRQ(ierr); 192 ierr = PetscObjectCompose(obj,"GetArrayPtr",0);CHKERRQ(ierr); 193 } else { 194 *lx = base + addr; 195 } 196 return 0; 197 } 198 199 #if defined(PETSC_HAVE_FORTRAN_CAPS) 200 #define petscisinfornanscalar_ PETSCISINFORNANSCALAR 201 #define petscisinfornanreal_ PETSCISINFORNANREAL 202 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 203 #define petscisinfornanscalar_ petscisinfornanscalar 204 #define petscisinfornanreal_ petscisinfornanreal 205 #endif 206 207 PETSC_EXTERN PetscBool petscisinfornanscalar_(PetscScalar *v) 208 { 209 return (PetscBool) PetscIsInfOrNanScalar(*v); 210 } 211 212 PETSC_EXTERN PetscBool petscisinfornanreal_(PetscReal *v) 213 { 214 return (PetscBool) PetscIsInfOrNanReal(*v); 215 } 216 217 218 219