GNU Octave/Octファイル

Octファイル

Octファイル(Octfile)は、GNU OctaveからC++の実行ファイルを呼び出す仕組みです。 Oct filesを使えば、自分で作ったC++の関数を、GNU Octaveからライブラリ関数とまったく同様に呼び出して使うことができます。 GNU Octaveに必要とする関数がない場合や、GNU Octaveで作った関数が実行速度や計算精度で物足りない場合に有用です。以下にOctfileのsampleを示します。

Octファイルの実行例

Hello world!

#include <octave/oct.h>

DEFUN_DLD (helloworld, args, nargout,
           "Hello World Help String")
{
  int nargin = args.length ();

  octave_stdout << "Hello World has "
                << nargin << " input arguments and "
                << nargout << " output arguments.\n";

  return octave_value_list ();
}

実数入力可算

//oct_AddTwoValues.cc
#include <octave/oct.h>
DEFUN_DLD(oct_AddTwoValues, args, ,
            "z=oct_AddTwoValues(x,y);")
{
    ColumnVector x(args(0).vector_value());
    ColumnVector y(args(1).vector_value());
    ColumnVector z;
    
    z = x + y;
    
    return octave_value (z); 
}

実行例

上記のC++プログラムを、DEFUN_DLDの最初の文字をファイル名として保存する。この例では、　oct_AddTwoValues.cc　として保存する。 GNU Octave上で以下を実行する。以下同様。

 GNU Octave, version 3.8.1
 Octave was configured for "i686-w64-mingw32".
>mkoctfile oct_AddTwoValues.cc
>oct_AddTwoValues(1,2)
⇒
ans =  3

A=(1:2)'
⇒
A =

   1
   2

B=(2:3)'
⇒
B =

   2
   3

oct_AddTwoValues(A,B)
⇒
   3
   5

ベクトル

//oct_ReturnTwoArrays.cc
#include <octave/oct.h>
DEFUN_DLD (oct_ReturnTwoArrays, args, ,
           "[C,D]=oct_ReturnTwoArrays(A,B);,A,B,C and D are arrays.\n")
{
  octave_value_list retval;
  int nargin = args.length ();
  if (nargin != 2) {
    print_usage();
  } else {
      //NDArray A = args(0).array_value();
      //NDArray B = args(1).array_value();
      NDArray A (args(0).array_value());
      NDArray B (args(1).array_value());
      if ( !error_state ) {
          retval(0) = A+B;
          retval(1) = A-B;
      }
  }
  return retval;
}

行列

//MatVecProd.cc
#include <octave/oct.h>
DEFUN_DLD (MatVecProd, args, , 
   "[nel,M,N,is_empty,ndims,Y]=MatVecProd(A,X); A is a matrix, X,Y are vectors.\n")
{
    octave_value_list retval;
    
    int nargin = args.length();
    if (nargin != 2){
        print_usage();
        return retval;
    }
    
    //Matrix A                 ( args(0).matrix_value());
    //Matrix X                 ( args(1).array_value());
    Matrix A                 = args(0).matrix_value();
    Matrix X                 = args(1).array_value();
    octave_idx_type nel      = A.numel();
    octave_idx_type M        = A.rows();    //Rows M
    octave_idx_type N        = A.columns(); //Colums N
    octave_idx_type is_empty = A.is_empty();
    octave_idx_type ndims    = A.ndims();
    Matrix Y(M,1);          //NDArray Y(M);　It didn't work!
    if ( error_state ) return retval;
    
    for ( int m=0; m < M; ++m ) {
        Y(m)=0;
        for ( int n=0; n < N; ++n ) {
            Y(m) += A(m,n) * X(n);
        }
    }
    
    retval(0)=nel;
    retval(1)=M;
    retval(2)=N;
    retval(3)=is_empty;
    retval(4)=ndims;
    retval(5)=Y;
    return retval;
}

GNU Octaveテストプログラム

%TEST_MatVecProd.m
function TEST_MatVecProd
clear all;
mkoctfile MatVecProd.cc
M=4
N=2
A=1:M*N;
A=reshape(A,M,N)
X=(1:N)'
[nel,M,N,is_empty,ndims,Y]=MatVecProd(A,X)
sizeY=size(Y)
AX=A*X
ERR=Y-AX

実行結果

　TEST_MatVecProd
⇒
M =  4
N =  2
A =
   1   5
   2   6
   3   7
   4   8
X =
   1
   2

nel =  8
M =  4
N =  2
is_empty = 0
ndims =  2
Y =
   11
   14
   17
   20

sizeY =
   4   1
AX =
   11
   14
   17
   20
ERR =
   0
   0
   0
   0

複素数

• Ref

//pippo.cc
#include <octave/oct.h>
DEFUN_DLD (pippo, args, , "test return statement")
{
    dim_vector dims;
    dims.resize (3);
    dims(2) = 4;
    dims(1) = 3;
    dims(0) = 2;

    ComplexNDArray a(dims, std::complex<double> (0.0, 0.0));
    // or 
    //Array<std::complex<double> > a (dims, std::complex<double> (0.0, 0.0));
    for (octave_idx_type i = 0; i < dims(0); i++) {
        for (octave_idx_type j = 0; j < dims(1); j++) {
            for (octave_idx_type k = 0; k < dims(2); k++) {
                a(i, j, k) = std::complex<double> (i, j + k);
            }
        }
    }
    return octave_value (a);
}

複素数入出力

//oct_complex.cc
#include <octave/oct.h>
DEFUN_DLD (oct_complex, args, , "[cx,cy,cz,Nc,Nr]=oct_complex(ca)\n")
{
    octave_value_list retval;
    int nargin = args.length();
    if (nargin != 1) {
        print_usage ();
        return retval;
    }
    
    ComplexMatrix cx(args(0).complex_matrix_value());//Input 0
    ComplexMatrix cy(args(0).complex_matrix_value());//for Output 1
    ComplexMatrix cz(args(0).complex_matrix_value());//for Output 2
    octave_idx_type  Nc = cx.columns();
    octave_idx_type  Nr = cx.rows();
    
    cy=real(cx);
    for (octave_idx_type i = 0; i < Nc; i++) {
        cz(i) = Complex(real(cx(i))*real(cx(i)),imag(cx(i))*imag(cx(i)));
    }
    
    retval(0)=cx;
    retval(1)=cy;
    retval(2)=cz;
    retval(3)=Nc;
    retval(4)=Nr;
    return retval;
}

実行例

 > mkoctfile oct_complex.cc
 > ca=[1+j 2+j];
 > [cx,cy,cz,Nc,Nr]=oct_complex(ca)
⇒ 
 cx =
    1 + 2i   3 + 4i
 cy =
    1   3
 cz =
    1 +  4i    9 + 16i
 Nc =  2
 Nr =  1

文字列

構造体

• RefA.1.5 Structures in Oct-Files

//oct_structdemo.cc
#include <octave/oct.h>
#include <octave/ov-struct.h>
DEFUN_DLD (oct_structdemo, args, , "Struct Demo")
{
    octave_value retval;
    int nargin = args.length ();
    if (args.length () == 2) {
        octave_scalar_map arg0 = args(0).scalar_map_value ();
        //octave_map arg0 = args(0).map_value();
        if ( ! error_state ) {
            std::string arg1 = args(1).string_value ();
            if ( ! error_state ) {
                octave_value tmp = arg0.contents (arg1);
                //octave_value tmp = arg0.contents (arg1)(0);
                if ( tmp.is_defined() ) {
                    octave_scalar_map st;
                    st.assign ("selected", tmp);
                    retval = octave_value (st);
                } else {
                    error("structdemo: struct does not have a field named '%s'\n",
                            arg1.c_str ());
                }
            } else {
                error ("structdemo: ARG2 must be a character string");
            }
        } else {
            error ("structdemo: ARG1 must be a struct");
        }
    } else {
        print_usage ();
    }
    return retval;
}

実行例

x.a = 1;
x.b="test";
x.c=[1, 2];
oct_structdemo (x, "b")
⇒ selected = test

セル配列

• Ref A.1.4 Cell Arrays in Oct-Files

//oct_celldemo.cc
#include <octave/oct.h>
#include <octave/Cell.h>
DEFUN_DLD (oct_celldemo, args, ,
            "b=oct_celldemo(c);,c is a cell, b are elements.")
{
    octave_value_list retval;
    int nargin = args.length();
    if (nargin != 1) {
        print_usage();
    } else {
        Cell c = args(0).cell_value();
        if ( ! error_state ) {
            for ( octave_idx_type i = 0; i < c.numel(); ++i ) {
                retval(i) = c(i);          // using operator syntax
                //retval(i) = c.elem (i);  // using method syntax
            }
        }
    }
    return retval;
}

実行例

[b1, b2, b3] = oct_celldemo({1, [1, 2], "test"})
⇒
b1 =  1
b2 =

   1   2

b3 = test

---

Copyright (C) 1996-2014 John W. Eaton.
Copyright (C) 2005-2007 MASUDA, Yutaka (日本語訳)
Copyright (C) 2014 Bethlehem4.
Permission is granted to make and distribute verbatim copies of this manual provided the copyright notice and this permission notice are preserved on all copies. Permission is granted to copy and distribute modified versions of this manual under the conditions for verbatim copying, provided that the entire resulting derived work is distributed under the terms of a permission notice identical to this one. Permission is granted to copy and distribute translations of this manual into another language, under the above conditions for modified versions.