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 F90Array1dCreateFortranAddr(array,start,len1,ptr)
36: implicit none
37: PetscInt start,len1
38: PetscFortranAddr, target :: &
39: & array(start:start+len1-1)
40: PetscFortranAddr, pointer :: ptr(:)
42: ptr => array
43: end subroutine
45: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
46: subroutine F90Array1dAccessScalar(ptr,address)
47: implicit none
48: PetscScalar, pointer :: ptr(:)
49: PetscFortranAddr address
50: PetscInt start
52: if (associated(ptr) .eqv. .false.) then
53: address = 0
54: else
55: start = lbound(ptr,1)
56: call F90Array1dGetAddrScalar(ptr(start),address)
57: endif
58: end subroutine
60: subroutine F90Array1dAccessReal(ptr,address)
61: implicit none
62: PetscReal, pointer :: ptr(:)
63: PetscFortranAddr address
64: PetscInt start
66: if (associated(ptr) .eqv. .false.) then
67: address = 0
68: else
69: start = lbound(ptr,1)
70: call F90Array1dGetAddrReal(ptr(start),address)
71: endif
72: end subroutine
74: subroutine F90Array1dAccessInt(ptr,address)
75: implicit none
76: PetscInt, pointer :: ptr(:)
77: PetscFortranAddr address
78: PetscInt start
80: if (associated(ptr) .eqv. .false.) then
81: address = 0
82: else
83: start = lbound(ptr,1)
84: call F90Array1dGetAddrInt(ptr(start),address)
85: endif
86: end subroutine
88: subroutine F90Array1dAccessFortranAddr(ptr,address)
89: implicit none
90: PetscFortranAddr, pointer :: ptr(:)
91: PetscFortranAddr address
92: PetscInt start
94: if (associated(ptr) .eqv. .false.) then
95: address = 0
96: else
97: start = lbound(ptr,1)
98: call F90Array1dGetAddrFortranAddr(ptr(start),address)
99: endif
100: end subroutine
102: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
103: subroutine F90Array1dDestroyScalar(ptr)
104: implicit none
105: PetscScalar, pointer :: ptr(:)
107: nullify(ptr)
108: end subroutine
110: subroutine F90Array1dDestroyReal(ptr)
111: implicit none
112: PetscReal, pointer :: ptr(:)
114: nullify(ptr)
115: end subroutine
117: subroutine F90Array1dDestroyInt(ptr)
118: implicit none
119: PetscInt, pointer :: ptr(:)
121: nullify(ptr)
122: end subroutine
124: subroutine F90Array1dDestroyFortranAddr(ptr)
125: implicit none
126: PetscFortranAddr, pointer :: ptr(:)
128: nullify(ptr)
129: end subroutine
130: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
131: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
132: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
133: subroutine F90Array2dCreateScalar(array,start1,len1, &
134: & start2,len2,ptr)
135: implicit none
136: PetscInt start1,len1
137: PetscInt start2,len2
138: PetscScalar, target :: &
139: & array(start1:start1+len1-1,start2:start2+len2-1)
140: PetscScalar, pointer :: ptr(:,:)
142: ptr => array
143: end subroutine
145: subroutine F90Array2dCreateReal(array,start1,len1, &
146: & start2,len2,ptr)
147: implicit none
148: PetscInt start1,len1
149: PetscInt start2,len2
150: PetscReal, target :: &
151: & array(start1:start1+len1-1,start2:start2+len2-1)
152: PetscReal, pointer :: ptr(:,:)
154: ptr => array
155: end subroutine
157: subroutine F90Array2dCreateInt(array,start1,len1, &
158: & start2,len2,ptr)
159: implicit none
160: PetscInt start1,len1
161: PetscInt start2,len2
162: PetscInt, target :: &
163: & array(start1:start1+len1-1,start2:start2+len2-1)
164: PetscInt, pointer :: ptr(:,:)
166: ptr => array
167: end subroutine
169: subroutine F90Array2dCreateFortranAddr(array,start1,len1, &
170: & start2,len2,ptr)
171: implicit none
172: PetscInt start1,len1
173: PetscInt start2,len2
174: PetscFortranAddr, target :: &
175: & array(start1:start1+len1-1,start2:start2+len2-1)
176: PetscFortranAddr, pointer :: ptr(:,:)
178: ptr => array
179: end subroutine
181: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
182: subroutine F90Array2dAccessScalar(ptr,address)
183: implicit none
184: PetscScalar, pointer :: ptr(:,:)
185: PetscFortranAddr address
186: PetscInt start1,start2
188: start1 = lbound(ptr,1)
189: start2 = lbound(ptr,2)
190: call F90Array2dGetAddrScalar(ptr(start1,start2),address)
191: end subroutine
193: subroutine F90Array2dAccessReal(ptr,address)
194: implicit none
195: PetscReal, pointer :: ptr(:,:)
196: PetscFortranAddr address
197: PetscInt start1,start2
199: start1 = lbound(ptr,1)
200: start2 = lbound(ptr,2)
201: call F90Array2dGetAddrReal(ptr(start1,start2),address)
202: end subroutine
204: subroutine F90Array2dAccessInt(ptr,address)
205: implicit none
206: PetscInt, pointer :: ptr(:,:)
207: PetscFortranAddr address
208: PetscInt start1,start2
210: start1 = lbound(ptr,1)
211: start2 = lbound(ptr,2)
212: call F90Array2dGetAddrInt(ptr(start1,start2),address)
213: end subroutine
215: subroutine F90Array2dAccessFortranAddr(ptr,address)
216: implicit none
217: PetscFortranAddr, pointer :: ptr(:,:)
218: PetscFortranAddr address
219: PetscInt start1,start2
221: start1 = lbound(ptr,1)
222: start2 = lbound(ptr,2)
223: call F90Array2dGetAddrFortranAddr(ptr(start1,start2),address)
224: end subroutine
226: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
227: subroutine F90Array2dDestroyScalar(ptr)
228: implicit none
229: PetscScalar, pointer :: ptr(:,:)
231: nullify(ptr)
232: end subroutine
234: subroutine F90Array2dDestroyReal(ptr)
235: implicit none
236: PetscReal, pointer :: ptr(:,:)
238: nullify(ptr)
239: end subroutine
241: subroutine F90Array2dDestroyInt(ptr)
242: implicit none
243: PetscInt, pointer :: ptr(:,:)
245: nullify(ptr)
246: end subroutine
248: subroutine F90Array2dDestroyFortranAddr(ptr)
249: implicit none
250: PetscFortranAddr, pointer :: ptr(:,:)
252: nullify(ptr)
253: end subroutine
254: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
255: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
256: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
257: subroutine F90Array3dCreateScalar(array,start1,len1, &
258: & start2,len2,start3,len3,ptr)
259: implicit none
260: PetscInt start1,len1
261: PetscInt start2,len2
262: PetscInt start3,len3
263: PetscScalar, target :: &
264: & array(start1:start1+len1-1,start2:start2+len2-1, &
265: & start3:start3+len3-1)
266: PetscScalar, pointer :: ptr(:,:,:)
268: ptr => array
269: end subroutine
271: subroutine F90Array3dCreateReal(array,start1,len1, &
272: & start2,len2,start3,len3,ptr)
273: implicit none
274: PetscInt start1,len1
275: PetscInt start2,len2
276: PetscInt start3,len3
277: PetscReal, target :: &
278: & array(start1:start1+len1-1,start2:start2+len2-1, &
279: & start3:start3+len3-1)
280: PetscReal, pointer :: ptr(:,:,:)
282: ptr => array
283: end subroutine
285: subroutine F90Array3dCreateInt(array,start1,len1, &
286: & start2,len2,start3,len3,ptr)
287: implicit none
288: PetscInt start1,len1
289: PetscInt start2,len2
290: PetscInt start3,len3
291: PetscInt, target :: &
292: & array(start1:start1+len1-1,start2:start2+len2-1, &
293: & start3:start3+len3-1)
294: PetscInt, pointer :: ptr(:,:,:)
296: ptr => array
297: end subroutine
299: subroutine F90Array3dCreateFortranAddr(array,start1,len1, &
300: & start2,len2,start3,len3,ptr)
301: implicit none
302: PetscInt start1,len1
303: PetscInt start2,len2
304: PetscInt start3,len3
305: PetscFortranAddr, target :: &
306: & array(start1:start1+len1-1,start2:start2+len2-1, &
307: & start3:start3+len3-1)
308: PetscFortranAddr, pointer :: ptr(:,:,:)
310: ptr => array
311: end subroutine
313: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
314: subroutine F90Array3dAccessScalar(ptr,address)
315: implicit none
316: PetscScalar, pointer :: ptr(:,:,:)
317: PetscFortranAddr address
318: PetscInt start1,start2,start3
320: start1 = lbound(ptr,1)
321: start2 = lbound(ptr,2)
322: start3 = lbound(ptr,3)
323: call F90Array3dGetAddrScalar(ptr(start1,start2,start3),address)
324: end subroutine
326: subroutine F90Array3dAccessReal(ptr,address)
327: implicit none
328: PetscReal, pointer :: ptr(:,:,:)
329: PetscFortranAddr address
330: PetscInt start1,start2,start3
332: start1 = lbound(ptr,1)
333: start2 = lbound(ptr,2)
334: start3 = lbound(ptr,3)
335: call F90Array3dGetAddrReal(ptr(start1,start2,start3),address)
336: end subroutine
338: subroutine F90Array3dAccessInt(ptr,address)
339: implicit none
340: PetscInt, pointer :: ptr(:,:,:)
341: PetscFortranAddr address
342: PetscInt start1,start2,start3
344: start1 = lbound(ptr,1)
345: start2 = lbound(ptr,2)
346: start3 = lbound(ptr,3)
347: call F90Array3dGetAddrInt(ptr(start1,start2,start3),address)
348: end subroutine
350: subroutine F90Array3dAccessFortranAddr(ptr,address)
351: implicit none
352: PetscFortranAddr, pointer :: ptr(:,:,:)
353: PetscFortranAddr address
354: PetscInt start1,start2,start3
356: start1 = lbound(ptr,1)
357: start2 = lbound(ptr,2)
358: start3 = lbound(ptr,3)
359: call F90Array3dGetAddrFortranAddr(ptr(start1,start2,start3), &
360: & address)
361: end subroutine
363: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
364: subroutine F90Array3dDestroyScalar(ptr)
365: implicit none
366: PetscScalar, pointer :: ptr(:,:,:)
368: nullify(ptr)
369: end subroutine
371: subroutine F90Array3dDestroyReal(ptr)
372: implicit none
373: PetscReal, pointer :: ptr(:,:,:)
375: nullify(ptr)
376: end subroutine
378: subroutine F90Array3dDestroyInt(ptr)
379: implicit none
380: PetscInt, pointer :: ptr(:,:,:)
382: nullify(ptr)
383: end subroutine
385: subroutine F90Array3dDestroyFortranAddr(ptr)
386: implicit none
387: PetscFortranAddr, pointer :: ptr(:,:,:)
389: nullify(ptr)
390: end subroutine
392: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
393: subroutine F90Array4dCreateScalar(array,start1,len1, &
394: & start2,len2,start3,len3,start4,len4,ptr)
395: implicit none
396: PetscInt start1,len1
397: PetscInt start2,len2
398: PetscInt start3,len3
399: PetscInt start4,len4
400: PetscScalar, target :: &
401: & array(start1:start1+len1-1,start2:start2+len2-1, &
402: & start3:start3+len3-1,start4:start4+len4-1)
403: PetscScalar, pointer :: ptr(:,:,:,:)
405: ptr => array
406: end subroutine
408: subroutine F90Array4dCreateReal(array,start1,len1, &
409: & start2,len2,start3,len3,start4,len4,ptr)
410: implicit none
411: PetscInt start1,len1
412: PetscInt start2,len2
413: PetscInt start3,len3
414: PetscInt start4,len4
415: PetscReal, target :: &
416: & array(start1:start1+len1-1,start2:start2+len2-1, &
417: & start3:start3+len3-1,start4:start4+len4-1)
418: PetscReal, pointer :: ptr(:,:,:,:)
420: ptr => array
421: end subroutine
423: subroutine F90Array4dCreateInt(array,start1,len1, &
424: & start2,len2,start3,len3,start4,len4,ptr)
425: implicit none
426: PetscInt start1,len1
427: PetscInt start2,len2
428: PetscInt start3,len3
429: PetscInt start4,len4
430: PetscInt, target :: &
431: & array(start1:start1+len1-1,start2:start2+len2-1, &
432: & start3:start3+len3-1,start4:start4+len4-1)
433: PetscInt, pointer :: ptr(:,:,:,:)
435: ptr => array
436: end subroutine
438: subroutine F90Array4dCreateFortranAddr(array,start1,len1, &
439: & start2,len2,start3,len3,start4,len4,ptr)
440: implicit none
441: PetscInt start1,len1
442: PetscInt start2,len2
443: PetscInt start3,len3
444: PetscInt start4,len4
445: PetscFortranAddr, target :: &
446: & array(start1:start1+len1-1,start2:start2+len2-1, &
447: & start3:start3+len3-1,start4:start4+len4-1)
448: PetscFortranAddr, pointer :: ptr(:,:,:,:)
450: ptr => array
451: end subroutine
453: subroutine F90Array4dAccessScalar(ptr,address)
454: implicit none
455: PetscScalar, pointer :: ptr(:,:,:,:)
456: PetscFortranAddr address
457: PetscInt start1,start2,start3,start4
459: start1 = lbound(ptr,1)
460: start2 = lbound(ptr,2)
461: start3 = lbound(ptr,3)
462: start4 = lbound(ptr,4)
463: call F90Array4dGetAddrScalar(ptr(start1,start2,start3,start4), &
464: & address)
465: end subroutine
467: subroutine F90Array4dAccessReal(ptr,address)
468: implicit none
469: PetscReal, pointer :: ptr(:,:,:,:)
470: PetscFortranAddr address
471: PetscInt start1,start2,start3,start4
473: start1 = lbound(ptr,1)
474: start2 = lbound(ptr,2)
475: start3 = lbound(ptr,3)
476: start4 = lbound(ptr,4)
477: call F90Array4dGetAddrReal(ptr(start1,start2,start3,start4), &
478: & address)
479: end subroutine
481: subroutine F90Array4dAccessInt(ptr,address)
482: implicit none
483: PetscInt, pointer :: ptr(:,:,:,:)
484: PetscFortranAddr address
485: PetscInt start1,start2,start3,start4
487: start1 = lbound(ptr,1)
488: start2 = lbound(ptr,2)
489: start3 = lbound(ptr,3)
490: start4 = lbound(ptr,4)
491: call F90Array4dGetAddrInt(ptr(start1,start2,start3,start4), &
492: & address)
493: end subroutine
495: subroutine F90Array4dAccessFortranAddr(ptr,address)
496: implicit none
497: PetscScalar, pointer :: ptr(:,:,:,:)
498: PetscFortranAddr address
499: PetscFortranAddr start1,start2,start3,start4
501: start1 = lbound(ptr,1)
502: start2 = lbound(ptr,2)
503: start3 = lbound(ptr,3)
504: start4 = lbound(ptr,4)
505: call F90Array4dGetAddrFortranAddr(ptr(start1,start2,start3, &
506: & start4),address)
507: end subroutine
509: subroutine F90Array4dDestroyScalar(ptr)
510: implicit none
511: PetscScalar, pointer :: ptr(:,:,:,:)
513: nullify(ptr)
514: end subroutine
516: subroutine F90Array4dDestroyReal(ptr)
517: implicit none
518: PetscReal, pointer :: ptr(:,:,:,:)
520: nullify(ptr)
521: end subroutine
523: subroutine F90Array4dDestroyInt(ptr)
524: implicit none
525: PetscInt, pointer :: ptr(:,:,:,:)
527: nullify(ptr)
528: end subroutine
530: subroutine F90Array4dDestroyFortranAddr(ptr)
531: implicit none
532: PetscFortranAddr, pointer :: ptr(:,:,:,:)
534: nullify(ptr)
535: end subroutine
537: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
538: !XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX!
539: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!