Actual source code: f90_fwrap.F90
1: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
2: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
3: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
4: #include <petsc/finclude/petscsys.h>
5: subroutine F90Array1dCreateScalar(array,start,len1,ptr)
6: implicit none
7: PetscInt start,len1
8: PetscScalar, target :: &
9: & array(start:start+len1-1)
10: PetscScalar, pointer :: ptr(:)
12: ptr => array
13: end subroutine
15: subroutine F90Array1dCreateReal(array,start,len1,ptr)
16: implicit none
17: PetscInt start,len1
18: PetscReal, target :: &
19: & array(start:start+len1-1)
20: PetscReal, pointer :: ptr(:)
22: ptr => array
23: end subroutine
25: subroutine F90Array1dCreateInt(array,start,len1,ptr)
26: implicit none
27: PetscInt start,len1
28: PetscInt, target :: &
29: & array(start:start+len1-1)
30: PetscInt, pointer :: ptr(:)
32: ptr => array
33: end subroutine
35: subroutine F90Array1dCreateMPIInt(array,start,len1,ptr)
36: implicit none
37: PetscInt start,len1
38: PetscMPIInt, target :: &
39: & array(start:start+len1-1)
40: PetscMPIInt, pointer :: ptr(:)
42: ptr => array
43: end subroutine
45: subroutine F90Array1dCreateFortranAddr(array,start,len1,ptr)
46: implicit none
47: PetscInt start,len1
48: PetscFortranAddr, target :: &
49: & array(start:start+len1-1)
50: PetscFortranAddr, pointer :: ptr(:)
52: ptr => array
53: end subroutine
55: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
56: subroutine F90Array1dAccessScalar(ptr,address)
57: implicit none
58: PetscScalar, pointer :: ptr(:)
59: PetscFortranAddr address
60: PetscInt start
62: if (associated(ptr) .eqv. .false.) then
63: address = 0
64: else
65: start = lbound(ptr,1)
66: call F90Array1dGetAddrScalar(ptr(start),address)
67: endif
68: end subroutine
70: subroutine F90Array1dAccessReal(ptr,address)
71: implicit none
72: PetscReal, pointer :: ptr(:)
73: PetscFortranAddr address
74: PetscInt start
76: if (associated(ptr) .eqv. .false.) then
77: address = 0
78: else
79: start = lbound(ptr,1)
80: call F90Array1dGetAddrReal(ptr(start),address)
81: endif
82: end subroutine
84: subroutine F90Array1dAccessInt(ptr,address)
85: implicit none
86: PetscInt, pointer :: ptr(:)
87: PetscFortranAddr address
88: PetscInt start
90: if (associated(ptr) .eqv. .false.) then
91: address = 0
92: else
93: start = lbound(ptr,1)
94: call F90Array1dGetAddrInt(ptr(start),address)
95: endif
96: end subroutine
98: subroutine F90Array1dAccessMPIInt(ptr,address)
99: implicit none
100: PetscMPIInt, pointer :: ptr(:)
101: PetscFortranAddr address
102: PetscInt start
104: if (associated(ptr) .eqv. .false.) then
105: address = 0
106: else
107: start = lbound(ptr,1)
108: call F90Array1dGetAddrMPIInt(ptr(start),address)
109: endif
110: end subroutine
112: subroutine F90Array1dAccessFortranAddr(ptr,address)
113: implicit none
114: PetscFortranAddr, pointer :: ptr(:)
115: PetscFortranAddr address
116: PetscInt start
118: if (associated(ptr) .eqv. .false.) then
119: address = 0
120: else
121: start = lbound(ptr,1)
122: call F90Array1dGetAddrFortranAddr(ptr(start),address)
123: endif
124: end subroutine
126: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
127: subroutine F90Array1dDestroyScalar(ptr)
128: implicit none
129: PetscScalar, pointer :: ptr(:)
131: nullify(ptr)
132: end subroutine
134: subroutine F90Array1dDestroyReal(ptr)
135: implicit none
136: PetscReal, pointer :: ptr(:)
138: nullify(ptr)
139: end subroutine
141: subroutine F90Array1dDestroyInt(ptr)
142: implicit none
143: PetscInt, pointer :: ptr(:)
145: nullify(ptr)
146: end subroutine
148: subroutine F90Array1dDestroyMPIInt(ptr)
149: implicit none
150: PetscMPIInt, pointer :: ptr(:)
152: nullify(ptr)
153: end subroutine
155: subroutine F90Array1dDestroyFortranAddr(ptr)
156: implicit none
157: PetscFortranAddr, pointer :: ptr(:)
159: nullify(ptr)
160: end subroutine
161: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
162: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
163: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
164: subroutine F90Array2dCreateScalar(array,start1,len1, &
165: & start2,len2,ptr)
166: implicit none
167: PetscInt start1,len1
168: PetscInt start2,len2
169: PetscScalar, target :: &
170: & array(start1:start1+len1-1,start2:start2+len2-1)
171: PetscScalar, pointer :: ptr(:,:)
173: ptr => array
174: end subroutine
176: subroutine F90Array2dCreateReal(array,start1,len1, &
177: & start2,len2,ptr)
178: implicit none
179: PetscInt start1,len1
180: PetscInt start2,len2
181: PetscReal, target :: &
182: & array(start1:start1+len1-1,start2:start2+len2-1)
183: PetscReal, pointer :: ptr(:,:)
185: ptr => array
186: end subroutine
188: subroutine F90Array2dCreateInt(array,start1,len1, &
189: & start2,len2,ptr)
190: implicit none
191: PetscInt start1,len1
192: PetscInt start2,len2
193: PetscInt, target :: &
194: & array(start1:start1+len1-1,start2:start2+len2-1)
195: PetscInt, pointer :: ptr(:,:)
197: ptr => array
198: end subroutine
200: subroutine F90Array2dCreateFortranAddr(array,start1,len1, &
201: & start2,len2,ptr)
202: implicit none
203: PetscInt start1,len1
204: PetscInt start2,len2
205: PetscFortranAddr, target :: &
206: & array(start1:start1+len1-1,start2:start2+len2-1)
207: PetscFortranAddr, pointer :: ptr(:,:)
209: ptr => array
210: end subroutine
212: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
213: subroutine F90Array2dAccessScalar(ptr,address)
214: implicit none
215: PetscScalar, pointer :: ptr(:,:)
216: PetscFortranAddr address
217: PetscInt start1,start2
219: start1 = lbound(ptr,1)
220: start2 = lbound(ptr,2)
221: call F90Array2dGetAddrScalar(ptr(start1,start2),address)
222: end subroutine
224: subroutine F90Array2dAccessReal(ptr,address)
225: implicit none
226: PetscReal, pointer :: ptr(:,:)
227: PetscFortranAddr address
228: PetscInt start1,start2
230: start1 = lbound(ptr,1)
231: start2 = lbound(ptr,2)
232: call F90Array2dGetAddrReal(ptr(start1,start2),address)
233: end subroutine
235: subroutine F90Array2dAccessInt(ptr,address)
236: implicit none
237: PetscInt, pointer :: ptr(:,:)
238: PetscFortranAddr address
239: PetscInt start1,start2
241: start1 = lbound(ptr,1)
242: start2 = lbound(ptr,2)
243: call F90Array2dGetAddrInt(ptr(start1,start2),address)
244: end subroutine
246: subroutine F90Array2dAccessFortranAddr(ptr,address)
247: implicit none
248: PetscFortranAddr, pointer :: ptr(:,:)
249: PetscFortranAddr address
250: PetscInt start1,start2
252: start1 = lbound(ptr,1)
253: start2 = lbound(ptr,2)
254: call F90Array2dGetAddrFortranAddr(ptr(start1,start2),address)
255: end subroutine
257: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
258: subroutine F90Array2dDestroyScalar(ptr)
259: implicit none
260: PetscScalar, pointer :: ptr(:,:)
262: nullify(ptr)
263: end subroutine
265: subroutine F90Array2dDestroyReal(ptr)
266: implicit none
267: PetscReal, pointer :: ptr(:,:)
269: nullify(ptr)
270: end subroutine
272: subroutine F90Array2dDestroyInt(ptr)
273: implicit none
274: PetscInt, pointer :: ptr(:,:)
276: nullify(ptr)
277: end subroutine
279: subroutine F90Array2dDestroyFortranAddr(ptr)
280: implicit none
281: PetscFortranAddr, pointer :: ptr(:,:)
283: nullify(ptr)
284: end subroutine
285: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
286: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
287: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
288: subroutine F90Array3dCreateScalar(array,start1,len1, &
289: & start2,len2,start3,len3,ptr)
290: implicit none
291: PetscInt start1,len1
292: PetscInt start2,len2
293: PetscInt start3,len3
294: PetscScalar, target :: &
295: & array(start1:start1+len1-1,start2:start2+len2-1, &
296: & start3:start3+len3-1)
297: PetscScalar, pointer :: ptr(:,:,:)
299: ptr => array
300: end subroutine
302: subroutine F90Array3dCreateReal(array,start1,len1, &
303: & start2,len2,start3,len3,ptr)
304: implicit none
305: PetscInt start1,len1
306: PetscInt start2,len2
307: PetscInt start3,len3
308: PetscReal, target :: &
309: & array(start1:start1+len1-1,start2:start2+len2-1, &
310: & start3:start3+len3-1)
311: PetscReal, pointer :: ptr(:,:,:)
313: ptr => array
314: end subroutine
316: subroutine F90Array3dCreateInt(array,start1,len1, &
317: & start2,len2,start3,len3,ptr)
318: implicit none
319: PetscInt start1,len1
320: PetscInt start2,len2
321: PetscInt start3,len3
322: PetscInt, target :: &
323: & array(start1:start1+len1-1,start2:start2+len2-1, &
324: & start3:start3+len3-1)
325: PetscInt, pointer :: ptr(:,:,:)
327: ptr => array
328: end subroutine
330: subroutine F90Array3dCreateFortranAddr(array,start1,len1, &
331: & start2,len2,start3,len3,ptr)
332: implicit none
333: PetscInt start1,len1
334: PetscInt start2,len2
335: PetscInt start3,len3
336: PetscFortranAddr, target :: &
337: & array(start1:start1+len1-1,start2:start2+len2-1, &
338: & start3:start3+len3-1)
339: PetscFortranAddr, pointer :: ptr(:,:,:)
341: ptr => array
342: end subroutine
344: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
345: subroutine F90Array3dAccessScalar(ptr,address)
346: implicit none
347: PetscScalar, pointer :: ptr(:,:,:)
348: PetscFortranAddr address
349: PetscInt start1,start2,start3
351: start1 = lbound(ptr,1)
352: start2 = lbound(ptr,2)
353: start3 = lbound(ptr,3)
354: call F90Array3dGetAddrScalar(ptr(start1,start2,start3),address)
355: end subroutine
357: subroutine F90Array3dAccessReal(ptr,address)
358: implicit none
359: PetscReal, pointer :: ptr(:,:,:)
360: PetscFortranAddr address
361: PetscInt start1,start2,start3
363: start1 = lbound(ptr,1)
364: start2 = lbound(ptr,2)
365: start3 = lbound(ptr,3)
366: call F90Array3dGetAddrReal(ptr(start1,start2,start3),address)
367: end subroutine
369: subroutine F90Array3dAccessInt(ptr,address)
370: implicit none
371: PetscInt, pointer :: ptr(:,:,:)
372: PetscFortranAddr address
373: PetscInt start1,start2,start3
375: start1 = lbound(ptr,1)
376: start2 = lbound(ptr,2)
377: start3 = lbound(ptr,3)
378: call F90Array3dGetAddrInt(ptr(start1,start2,start3),address)
379: end subroutine
381: subroutine F90Array3dAccessFortranAddr(ptr,address)
382: implicit none
383: PetscFortranAddr, pointer :: ptr(:,:,:)
384: PetscFortranAddr address
385: PetscInt start1,start2,start3
387: start1 = lbound(ptr,1)
388: start2 = lbound(ptr,2)
389: start3 = lbound(ptr,3)
390: call F90Array3dGetAddrFortranAddr(ptr(start1,start2,start3), &
391: & address)
392: end subroutine
394: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
395: subroutine F90Array3dDestroyScalar(ptr)
396: implicit none
397: PetscScalar, pointer :: ptr(:,:,:)
399: nullify(ptr)
400: end subroutine
402: subroutine F90Array3dDestroyReal(ptr)
403: implicit none
404: PetscReal, pointer :: ptr(:,:,:)
406: nullify(ptr)
407: end subroutine
409: subroutine F90Array3dDestroyInt(ptr)
410: implicit none
411: PetscInt, pointer :: ptr(:,:,:)
413: nullify(ptr)
414: end subroutine
416: subroutine F90Array3dDestroyFortranAddr(ptr)
417: implicit none
418: PetscFortranAddr, pointer :: ptr(:,:,:)
420: nullify(ptr)
421: end subroutine
423: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
424: subroutine F90Array4dCreateScalar(array,start1,len1, &
425: & start2,len2,start3,len3,start4,len4,ptr)
426: implicit none
427: PetscInt start1,len1
428: PetscInt start2,len2
429: PetscInt start3,len3
430: PetscInt start4,len4
431: PetscScalar, target :: &
432: & array(start1:start1+len1-1,start2:start2+len2-1, &
433: & start3:start3+len3-1,start4:start4+len4-1)
434: PetscScalar, pointer :: ptr(:,:,:,:)
436: ptr => array
437: end subroutine
439: subroutine F90Array4dCreateReal(array,start1,len1, &
440: & start2,len2,start3,len3,start4,len4,ptr)
441: implicit none
442: PetscInt start1,len1
443: PetscInt start2,len2
444: PetscInt start3,len3
445: PetscInt start4,len4
446: PetscReal, target :: &
447: & array(start1:start1+len1-1,start2:start2+len2-1, &
448: & start3:start3+len3-1,start4:start4+len4-1)
449: PetscReal, pointer :: ptr(:,:,:,:)
451: ptr => array
452: end subroutine
454: subroutine F90Array4dCreateInt(array,start1,len1, &
455: & start2,len2,start3,len3,start4,len4,ptr)
456: implicit none
457: PetscInt start1,len1
458: PetscInt start2,len2
459: PetscInt start3,len3
460: PetscInt start4,len4
461: PetscInt, target :: &
462: & array(start1:start1+len1-1,start2:start2+len2-1, &
463: & start3:start3+len3-1,start4:start4+len4-1)
464: PetscInt, pointer :: ptr(:,:,:,:)
466: ptr => array
467: end subroutine
469: subroutine F90Array4dCreateFortranAddr(array,start1,len1, &
470: & start2,len2,start3,len3,start4,len4,ptr)
471: implicit none
472: PetscInt start1,len1
473: PetscInt start2,len2
474: PetscInt start3,len3
475: PetscInt start4,len4
476: PetscFortranAddr, target :: &
477: & array(start1:start1+len1-1,start2:start2+len2-1, &
478: & start3:start3+len3-1,start4:start4+len4-1)
479: PetscFortranAddr, pointer :: ptr(:,:,:,:)
481: ptr => array
482: end subroutine
484: subroutine F90Array4dAccessScalar(ptr,address)
485: implicit none
486: PetscScalar, pointer :: ptr(:,:,:,:)
487: PetscFortranAddr address
488: PetscInt start1,start2,start3,start4
490: start1 = lbound(ptr,1)
491: start2 = lbound(ptr,2)
492: start3 = lbound(ptr,3)
493: start4 = lbound(ptr,4)
494: call F90Array4dGetAddrScalar(ptr(start1,start2,start3,start4), &
495: & address)
496: end subroutine
498: subroutine F90Array4dAccessReal(ptr,address)
499: implicit none
500: PetscReal, pointer :: ptr(:,:,:,:)
501: PetscFortranAddr address
502: PetscInt start1,start2,start3,start4
504: start1 = lbound(ptr,1)
505: start2 = lbound(ptr,2)
506: start3 = lbound(ptr,3)
507: start4 = lbound(ptr,4)
508: call F90Array4dGetAddrReal(ptr(start1,start2,start3,start4), &
509: & address)
510: end subroutine
512: subroutine F90Array4dAccessInt(ptr,address)
513: implicit none
514: PetscInt, pointer :: ptr(:,:,:,:)
515: PetscFortranAddr address
516: PetscInt start1,start2,start3,start4
518: start1 = lbound(ptr,1)
519: start2 = lbound(ptr,2)
520: start3 = lbound(ptr,3)
521: start4 = lbound(ptr,4)
522: call F90Array4dGetAddrInt(ptr(start1,start2,start3,start4), &
523: & address)
524: end subroutine
526: subroutine F90Array4dAccessFortranAddr(ptr,address)
527: implicit none
528: PetscScalar, pointer :: ptr(:,:,:,:)
529: PetscFortranAddr address
530: PetscFortranAddr start1,start2,start3,start4
532: start1 = lbound(ptr,1)
533: start2 = lbound(ptr,2)
534: start3 = lbound(ptr,3)
535: start4 = lbound(ptr,4)
536: call F90Array4dGetAddrFortranAddr(ptr(start1,start2,start3, &
537: & start4),address)
538: end subroutine
540: subroutine F90Array4dDestroyScalar(ptr)
541: implicit none
542: PetscScalar, pointer :: ptr(:,:,:,:)
544: nullify(ptr)
545: end subroutine
547: subroutine F90Array4dDestroyReal(ptr)
548: implicit none
549: PetscReal, pointer :: ptr(:,:,:,:)
551: nullify(ptr)
552: end subroutine
554: subroutine F90Array4dDestroyInt(ptr)
555: implicit none
556: PetscInt, pointer :: ptr(:,:,:,:)
558: nullify(ptr)
559: end subroutine
561: subroutine F90Array4dDestroyFortranAddr(ptr)
562: implicit none
563: PetscFortranAddr, pointer :: ptr(:,:,:,:)
565: nullify(ptr)
566: end subroutine
568: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
569: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
570: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!