C++ interface to Fortran with pointers management

HySoP/CMakeLists.txt

@@ -97,6 +97,7 @@ endif(EXISTS ${CMAKE_SOURCE_DIR}/config.hpp.cmake)
 # The list of all dirs containing sources to be compiled
 set(${PROJECT_NAME}_SRCDIRS 
 src
+src/interfaces/Fortran2Cpp
 src/interfaces/ppm
 )
 # List of the possible file names extensions

HySoP/src/callPPM.f90

@@ -23,7 +23,29 @@ module testPPM
 
   implicit none
 
-  contains
+contains
+
+  subroutine mult(A,x,sizeA)
+    real, dimension(:) :: A
+    real :: x
+    integer :: sizeA,i,j
+    
+    real*8, dimension(sizeA) :: B
+    
+   print *, "otot"
+   print *, A(1)
+   print *, "titi"
+
+    do j=1,sizeA
+       B(j) =j
+       !! do i=1,sizeA
+   !!       print *, A(i,j)
+!!$          A(i,j) = x*A(i,j)
+       !!print *, A(j)
+       !!print *,'toto', B(j)
+       !!end do
+    end do
+  end subroutine mult
 
   subroutine init()
 

HySoP/src/interfaces/Fortran2Cpp/WrapC.hpp

+#ifndef CWRAPPER_HPP
+#define CWRAPPER_HPP 
+
+typedef struct {
+  int length;
+  double* elements;
+} C2FPtr;
+
+extern "C" {
+
+  void wrapC2F_allocatedPtr(double*, int*);
+  void wrapC2F_NULLPtr(C2FPtr*);
+  void wrapC2F_NULLPtrBis(double**, int*);
+ 
+}
+
+#endif

HySoP/src/interfaces/Fortran2Cpp/WrapFortran.f90

+module WrapFort
+
+  use  iso_c_binding
+  use modTest
+
+  implicit none
+
+  public wrapC2F_allocatedPtr, wrapC2F_NULLPtr, wrapC2F_NULLPtrBis 
+
+  !> A structure to bind C and fortran pointers 
+  type, bind(C) :: C2FPtr
+      integer (c_int) :: length
+      type (c_ptr) :: elements
+   end type C2FPtr
+
+   logical, parameter :: NDEBUG = .TRUE.
+
+contains
+
+  !> Send an already allocated C pointer to a Fortran subroutine.
+  !! @param type(c_Ptr) a C pointer (void*)
+  !! @param type(c_int) size of C pointer
+  subroutine wrapC2F_allocatedPtr(cptr, sizeCptr) bind(C, name='wrapC2F_allocatedPtr')
+
+    type(c_Ptr),intent(in),VALUE :: cptr
+    integer (kind=c_int), intent(IN) :: sizeCptr
+
+    real(kind=c_double), pointer, dimension(:) :: xp => NULL()
+
+    if(NDEBUG) print *, '=== wrapC2F_allocatedPtr ===' 
+    ! Associate cptr and xp. 
+    call c_f_pointer (cptr, xp, (/sizeCPtr/))
+
+    if(.not.associated(xp) ) then
+       print *, 'Error, association failed'
+    end if
+
+    ! Do some stuff on xp ...
+    !    xp(1) = -3.9
+    call cas3(xp)
+    print *, '=== End wrapC2F_allocatedPtr === '
+  end subroutine WrapC2F_allocatedPtr
+
+
+  !> Send a NULL C pointer to a Fortran subroutine and get it back properly allocated 
+  !! @param[inout] type(c_Ptr) a C pointer (void*)
+  !! @param[out] type(c_int) size of C pointer
+  subroutine wrapC2F_NULLPtr(vector) bind(C, name='wrapC2F_NULLPtr')
+
+    type(C2FPtr) :: vector
+
+    real(kind=c_double), pointer, dimension(:) :: xp => NULL()
+
+    if(NDEBUG) print *, '=== wrapC2F_NULLPtr ===' 
+
+    call cas5(xp)
+
+    if(.not.associated(xp) ) then
+       print *, 'Error, association failed'
+    end if
+
+    vector%length = size(xp)
+    vector%elements = c_loc(xp(1))
+
+    print *, '=== End of wrapC2F_NULLPtr === '
+  end subroutine wrapC2F_NULLPtr
+
+  !> Send a NULL C pointer to a Fortran subroutine and get it back properly allocated 
+  !! @param[inout] type(c_Ptr) a C pointer (void*)
+  !! @param[out] type(c_int) size of C pointer
+  subroutine wrapC2F_NULLPtrBis(vector, length) bind(C, name='wrapC2F_NULLPtrBis')
+
+    type(c_ptr),intent(inout) :: vector
+    integer(c_int), intent(out) :: length
+
+    real(kind=c_double), pointer, dimension(:) :: xp => NULL()
+
+    if(NDEBUG) print *, '=== wrapC2F_NULLPtr ===' 
+
+    call cas5(xp)
+
+!!$    if(.not.associated(xp) ) then
+!!$       print *, 'Error, association failed'
+!!$    end if
+!!$    length = size(xp)
+!!$    vector= c_loc(xp(1))
+
+    call aliasF2C(vector, xp, length)
+    
+    print *, '=== End of wrapC2F_NULLPtr === '
+  end subroutine wrapC2F_NULLPtrBis
+
+  subroutine aliasF2C(vectorC, vectorF, length)
+    type(c_ptr),intent(inout) :: vectorC
+    integer(c_int), intent(out) :: length
+    real(c_double), pointer, dimension(:) :: vectorF
+
+    if(.not.associated(vectorF) ) then
+       stop 'Error, input Fortran vector is not associated'
+    end if
+    length = size(vectorF)
+    vectorC = c_loc(vectorF(1))
+    
+  end subroutine aliasF2C
+  
+end module WrapFort

HySoP/src/interfaces/ppm/ppm_wrapper.hpp

@@ -5,24 +5,38 @@
 #include"FCMangle.h"
 #include<mpi.h>
 #include <cstring>
+#include "Mesh.hpp"
+
 /** Namespace for ppm functions and subroutines */
+
 namespace PPM 
 {
   
   extern "C" {
+
+    //PPM_GLOBAL(pass,PASS);
     // Init and close ppm
     void PPM_MODULE(ppm_module_init,ppm_init, PPM_MODULE_INIT,PPM_INIT)(int*, int*,int*,int*,int*,int*,int*,int*,int*);
     void PPM_MODULE(ppm_module_finalize, ppm_finalize,PPM_MODULE_FINALIZE,PPM_FINALIZE)(int&);
     // Display functions
     void PPM_MODULE(charfunctions, start, CHARFUNCTIONS, START)(double*,int*,char*,int );
     void PPM_MODULE(charfunctions, stop, CHARFUNCTIONS, stop)(double*,int*,char*,int );
+
+    void PPM_MODULE(testppm,mult,TESTPPM,MULT)(double*,double*,int*);
+
     // Topologies
-    //void PPM_MODULE(ppm_module_mktopo, ppm_mktopo, PPM_MODULE_MKTOPO,PPM_MKTOPO)(int*, 
+    //    void PPM_MODULE(ppm_module_mktopo, ppm_mktopo, PPM_MODULE_MKTOPO,PPM_MKTOPO)(int*, double*, int*,
     
+    void PPM_MODULE(modtest, cas1, MODTEST,CAS1)(double*);
+    void PPM_MODULE(modtest, cas2, MODTEST,CAS2)(double*);
+    void PPM_MODULE(modtest, cas3, MODTEST,CAS3)(double*);
+    void PPM_MODULE(modtest, cas4, MODTEST,CAS4)(double*);
+    void PPM_MODULE(modtest, cas5, MODTEST,CAS5)(double*);
+    void PPM_MODULE(modtest, cas6, MODTEST,CAS6)(double*);
+    void PPM_MODULE(modtest, application3, MODTEST,APPLICATION3)(); 
     //    int PPM_MODULE(ppm_module_data,ppm_debug, PPM_MODULE_DATA,PPM_DEBUG);
   }
 
-  //#define DEBUG PPM_MODULE_(ppm_module_data,ppm_debug, PPM_MODULE_DATA,PPM_DEBUG)
 
   /** A static class to call wrapped ppm functions */
   class wrapper 
@@ -86,4 +100,6 @@ namespace PPM
   };
 
 }
+
+
 #endif

HySoP/src/main.cxx

@@ -10,52 +10,152 @@
 #include<ParmesDef.hpp>
 #include<Domain.hpp>
 #include<mpi.h>
-using namespace std;
+#include "WrapC.hpp"
+
+using namespace std ;
 
 int main(int argc, char* argv[])
 {
-        float p1_x,p1_y,p2_x,p2_y;
+  float p1_x,p1_y,p2_x,p2_y;
 
-        p1_x = 0.0;
-        p1_y = 0.0;
-        p2_x = 3.0;
-        p2_y = 4.0;
+  p1_x = 0.0;
+  p1_y = 0.0;
+  p2_x = 3.0;
+  p2_y = 4.0;
 
-	double t0; 
-	int info;
+  double t0; 
+  int info;
+	
+  string str = "une chaîne de caractères";
+  // créer le buffer pour copier la chaîne
+  size_t size = str.size() + 1;
+  char * buffer = new char[ size ];
+  // copier la chaîne
+  strncpy( buffer, str.c_str(), size );
 	
-	string str = "une chaîne de caractères";
-	// créer le buffer pour copier la chaîne
-	size_t size = str.size() + 1;
-	char * buffer = new char[ size ];
-	// copier la chaîne
-	strncpy( buffer, str.c_str(), size );
+  //	__ppm_module_substart_MOD_substart(buffer, &t0, &info);
+
+
+  // =====  Physical domain definition =====
+  // dimensions 
+  Parmes::Def::vector3D dimsD = { { 1.0, 3.1, 4.3} };
+  // "Lowest" point
+  Parmes::Def::vector3D startPoint = { { 0., 1., 2.} };
+  // The domain
+  Parmes::Model::Domain<3> domain(dimsD, startPoint);
+  // =====  Grid definition =====
+  // Number of points in each dir ...
+  boost::array<size_t, 3>  nbSteps = { { 3, 4, 5} };
+  Parmes::Discr::Grid<3> grid(domain, nbSteps);
 	
-	//	__ppm_module_substart_MOD_substart(buffer, &t0, &info);
+  std::cout << grid << std::endl;
+
+  std::string msg = "Main Programm in c++";
+
+  MPI::Init();
+  MPI::Intracomm Comm = MPI::COMM_WORLD;
+  PPM::wrapper::init(3, 8, -15, Comm, 2, &info, 99, 98,97);
+  PPM::wrapper::substart(msg, &t0, &info);
+
+  MPI::Finalize();
+  PPM::wrapper::substop(msg, &t0, &info);
 
+  double A[1][1];
+  double x = 2.0;
+  double B[1];
+  int sizeA = 1;
+  for (int i = 0; i<sizeA; ++i)
+    for(int j =0; j<sizeA;++j)
+      A[i][j] = i+j;
+  B[0] = 12;
+  //B[1] = 4;
+  //B[2] =2;
 
-	// =====  Physical domain definition =====
-	// dimensions 
-	Parmes::Def::vector3D dimsD = { { 1.0, 3.1, 4.3} };
-	// "Lowest" point
-	Parmes::Def::vector3D startPoint = { { 0., 1., 2.} };
-	// The domain
-	Parmes::Model::Domain<3> domain(dimsD, startPoint);
-	// =====  Grid definition =====
-	// Number of points in each dir ...
-	boost::array<size_t, 3>  nbSteps = { { 3, 4, 5} };
-	Parmes::Discr::Grid<3> grid(domain, nbSteps);
 	
-	std::cout << grid << std::endl;
+ //  for (int i = 0; i<sizeA; ++i)
+//     {
+//       for(int j =0; j<sizeA;++j)
+// 	std::cout << A[i][j] << " ";
+//       std::cout << std::endl;
+//     }
+	
+  //	PPM::PPM_MODULE(testppm,mult,TESTPPM,MULT)(B,&x,&sizeA);
+
+//   for (int i = 0; i<sizeA; ++i)
+//     {
+//       std::cout<< "B " << B[i] << std::endl;
+//       for(int j =0; j<sizeA;++j)
+// 	std::cout << A[i][j] << " ";
+//       std::cout << std::endl;
+//     }
+
+  std::vector<double> vect(2);
+	
+  vect[0] = 12;
+  vect[1] =14.8;
+	
+  double * vect2 = NULL;
+  double * vect3;
+
+  int sizeV = 2;
+  //	PPM::PPM_MODULE(modtest, cas1, MODTEST,CAS1)(&vect[0]);
+  //PPM::PPM_MODULE(modtest, cas2, MODTEST,CAS2)(&vect[0]);//, &sizeV);
+  //PPM::PPM_MODULE(modtest, cas3, MODTEST,CAS3)(vect.data(), &sizeV);//, &sizeV);
+  //	PPM::PPM_MODULE(modtest, cas4, MODTEST,CAS4)(vect2);//, &sizeV);
+  //	PPM::PPM_MODULE(modtest, testglob, MODTEST,TESTGLOB)();
+  // 	std::cout << "iauziazuizauzi" << std::endl;
+  // 	double** adrVec;
+  // 	adrVec = new(double*);
+  // 	*adrVec = &vect[0];
+
+  // 	cout << " titi" << (*adrVec)[0] << endl;
+
+  //	PPM::PPM_MODULE(modtest, cas6, MODTEST,CAS6)(adrVec);//, &sizeV);
+  //	std::cout << "yiha" << vect2[1] << std::endl;
+	  
+  //Mesh* p_obj = 0 ;
+	
+  double * toto = 0;
+
+  int nx=2, ny=2, nz=1;
+  double toto2[nx];
+
+  toto2[0] = 12.8;
+  toto2[nx-1]= 90.8;
+  //PPM::wrapper::bridgeC2F(&nx,toto2);
+	
+  //PPM::PPM_MODULE(modtest, bridge, MODTEST,BRIDGE)(&nx,&ny,&nz,toto);
+  //PPM::PPM_MODULE(modtest, application3, MODTEST,APPLICATION3)();
+	
+  //PPM::PPM_MODULE(modtest, bridge2, MODTEST,BRIDGE2)(&nx,toto2);
+	
+  double * toto3 = 0; 
+  //PPM::PPM_MODULE(modtest, testglob, MODTEST,TESTGLOB)(toto,&toto3);
+  int ntoto;
+  //	PPM::PPM_MODULE(wrapfort, bridgef2c, WRAPFORTMODTEST,BRIDGEF2C)(&ntoto,toto);
+	
+
+  //getFortranPtr(&toto);
+	
+
+  //std::cout << "toto toto " << toto[0] << " " << toto[1]<< std::endl;
+  //	std::cout << "totototo " << toto2[0] << " " << toto2[1]<< std::endl;
 
-	std::string msg = "Main Programm in c++";
+  std::cout << "TEST PPM ..." << std::endl;
 
-	MPI::Init();
-	MPI::Intracomm Comm = MPI::COMM_WORLD;
-	PPM::wrapper::init(3, 8, -15, Comm, 2, &info, 99, 98,97);
-	PPM::wrapper::substart(msg, &t0, &info);
+  C2FPtr * myVector = new C2FPtr;
+  wrapC2F_NULLPtr(myVector);
+  
+  for(int i =0; i<myVector->length; ++i)
+    cout << myVector->elements[i] << endl;
+  
+  wrapC2F_NULLPtrBis(&toto, &ntoto);
+  
+  cout << ntoto << endl;
+  for(int i = 0; i<ntoto ;++i)
+    cout << toto[i] << endl;
+  
+  delete(myVector);
 
-	MPI::Finalize();
-	PPM::wrapper::substop(msg, &t0, &info);
 }
 

HySoP/src/modTest.f90

+module modTest
+
+  implicit none
+   
+contains
+
+  ! IN : an already allocated pointer with explicit size. 
+  subroutine cas1(x)
+
+    real*8, dimension(2) :: x
+
+    x(2) = x(2) +1.65
+    print *, 'cas 1', x(1), ' ', x(2)
+
+  end subroutine cas1
+
+  ! IN an already allocated pointer with implicit size. 
+  subroutine cas2(x)
+
+    real*8, dimension(*), intent(inout) :: x
+
+
+    print *, 'cas2a', x(1), ' ', x(2)
+
+    x(2) = x(2) +1.65
+    print *, 'cas 2b', x(1), ' ', x(2)
+
+  end subroutine cas2
+
+  ! Already allocated pointer, no size information, intent(IN)
+  subroutine cas3(x)
+
+    real(kind=8), dimension(:),  intent(in) :: x
+
+    !print *, 'cas3 ', shape(x) 
+    print *, 'cas3a ', x(1), ' ', x(2)
+
+  end subroutine cas3
+
+  subroutine cas4(x)
+
+    real(kind=8), dimension(:),  intent(inout) :: x
+
+
+    print *, 'cas3a', x(1), ' ', x(2)
+
+    x(2) = x(2) +1.65
+    print *, 'cas 3b', x(1), ' ', x(2)
+
+  end subroutine cas4
+
+  subroutine cas5(x)
+
+    !!   integer, intent(in) :: size
+    real(kind=8), pointer, dimension(:) :: x
+
+    allocate(x(2))
+    x(1) =12
+    x(2) = 8
+    print *, 'cas5a', x(1), ' ', x(2)
+
+    x(2) = 1.65
+    print *, 'cas 5b', x(1), ' ', x(2)
+    print * , x
+  end subroutine cas5
+
+  subroutine cas6(x)
+
+    !!   integer, intent(in) :: size
+    real(kind=8), pointer, dimension(:) :: x
+
+    x(1) =12
+    x(2) = 8
+    print *, 'cas6a', shape(x), ' ', x(1), ' ', x(2)
+
+    x(2) = 1.65
+    print *, 'cas 6b', x(1), ' ', x(2)
+    print * , x
+
+  end subroutine cas6
+
+  subroutine Application3()
+    use temp
+    real, dimension(:), pointer :: x
+    call bibi(x)
+    print *, x
+  end subroutine Application3
+
+end module modTest