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