Fortran

From Octave
Jump to navigation Jump to search
The printable version is no longer supported and may have rendering errors. Please update your browser bookmarks and please use the default browser print function instead.

This page describes an example of how to call liboctave functions from a Fortran program. In the example we will load a single matrix, stored in ASCII format, from a data file. It consists of two steps:

  1. write a C++ function with a C compatible interface and C linkage that reads a variable from an Octave ASCII file
  2. write Fortran code using the "iso_c_binding" intrinsic module to call the C++ function

Data file

File: data.txt
 1.59797350e-01 5.41307474e-01 1.12127655e-01 2.09249248e-01
 3.22564589e-01 7.94307947e-01 8.25924316e-01 5.38352076e-01
 3.63990736e-01 1.90371212e-02 2.89370865e-01 1.30131246e-01
 6.28360462e-01 1.98831330e-01 6.89539723e-01 6.91062597e-01

The file was generated with

A = rand (4);
save -ascii data.txt A

C++ function

C++ function to load a single matrix, stored in ASCII format, from a data file.

File: octave_file_io.cc
// Octave header
#include <octave/oct.h>
#include <octave/ls-mat-ascii.h>

// Custom header containing the C compatible interface
#include <octave_file_io.h>

//! Load a single matrix, stored in ASCII format, from a data file.
//!
//! @param file_name name of the data file.
//! @param data pointer to the read-in matrix stored as fortran vector
//!             (column-major order).
//! @param numel number of elements in @p data.

int octave_load (const char* file_name, double** data, int* numel)
{
  // Define variable to hold the read data.
  octave_value read_data;

  // Read a plain ASCII matrix from data file.
  std::ifstream in_file_stream (file_name, std::ios::binary);
  read_mat_ascii_data (in_file_stream, file_name, read_data);
  in_file_stream.close ();

  // Convert read data to numerical array (matrix).
  NDArray A = read_data.array_value ();

  // Extract number of elements in matrix A.
  *numel = A.numel ();

  // Allocate memory to pointer to returned values.
  *data  = (double*) malloc (A.numel () * sizeof (double));

  // Copy the content of matrix A to data structure Fortran can handle.
  memcpy (*data, A.fortran_vec (), A.numel () * sizeof (double));

  return 0;
}

Header file

Header file with C interface to octave_file_io.cc.

File: octave_file_io.h
#ifndef OCTAVE_FILE_IO_H
#define OCTAVE_FILE_IO_H

#ifdef __cplusplus
extern "C" {
#endif
  
  int octave_load (const char* filename, double** data, int* numel);

#ifdef __cplusplus
}
#endif

#endif

Fortran Code

Fortran main program to read the plain ASCII matrix with the help of the Octave-C++ function. The read in matrix is printed to the screen.

File: octave_file_io_example.f90
program octave_file_io_example

  use iso_c_binding

  implicit none

  interface
     function octave_load (filename, data, numel) bind(c, name="octave_load")

       use iso_c_binding
       implicit none

       integer(c_int) :: octave_load

       character(kind=c_char), intent(in) :: filename(*)

       type(c_ptr), intent(out) :: data
       integer(c_int), intent(out) :: numel

     end function octave_load
  end interface

  integer(c_int) :: res
  type(c_ptr) :: data
  real(c_double), pointer :: fdata(:)
  integer(c_int) :: numel

  res = octave_load (c_char_'data.txt' // c_null_char, data, numel)

  call c_f_pointer (data, fdata, (/numel/))

  write (*,*) fdata

end program octave_file_io_example

Compiling the code

Generate octave_file_io.o from octave_file_io.cc.

mkoctfile -I. -c octave_file_io.cc

Generate octave_file_io_example.exe from octave_file_io_example.f90 including octave_file_io.o.

mkoctfile -I. --link-stand-alone octave_file_io_example.f90 octave_file_io.o -o octave_file_io_example.exe -lgfortran

If you receive errors about missing libraries, make sure your LD_LIBRARY_PATH is set correctly to find all Octave libraries, e.g.

$ ldd ./octave_file_io_example.exe
...
       libgfortran.so.4 => /usr/lib64/libgfortran.so.4 (0x00007fe9eb62b000)
       liboctinterp.so.8 => not found
       liboctave.so.8 => not found
...

Then find liboctinterp.so.8 and liboctave.so.8 on your system and type

export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/path/to/lib/octave/9.1.0/