Re: compiling C++ and C functions for octave

[Iago Mosqueira]

OK. Many thanks again. It now works as wished.


Iago

----- Original Message -----
From: "Geraint" <address@hidden>
To: "Iago Mosqueira" <address@hidden>
Cc: <address@hidden>
Sent: 25 November 2002 18:43
Subject: Re: compiling C++ and C functions for octave


> "Iago Mosqueira" <address@hidden> writes:
>
> > Just one question. For this code to accept vectors, would it need to run
a
> > loop? Or can the octave C++ functions take care of that?
> >
> > Thanks,
> >
> > iago
> >
>
> It can accept vectors with a slight modification (you might want to
improve the error handling, check that all the arguments are the same
length, etc.). I don't know anything about the R mathlib, but if the R
functions accept vector inputs you could probably eliminate the for loop in
the DLD.
>
> /Geraint.
>
>
> octave:1> matern(1,2,3)
> ans = 0.96965
> octave:2> matern(2,3,4)
> ans = 0.96396
> octave:3> matern([1;2],[2;3],[3;4])
> ans =
>
>   0.96965
>   0.96396
>
>
> matern.cc: (compile with 'mkoctfile -lRmath matern.cc')
>
> #define MATHLIB_STANDALONE
> #include <Rmath.h>
> #include <math.h>
> #include <octave/oct.h>
>
> extern "C"
> {
>   float matern(float x, float phi, int kappa)
>   {
>     /*Returns the Matern function for x, phi and kappa.*/
>
>     /* Variables */
>     float ans, cte;
>     float uphi=x/phi;
>
>     /* Matern */
>
>     if (uphi==0) return 1;
>     else{
>       if (kappa==0.5)
> ans = exp(-uphi);
>       else {
> cte = R_pow(2, (-(kappa-1)))/gammafn(kappa);
> ans = cte * R_pow(uphi, kappa) * bessel_k(uphi,kappa,1);
>       }
>     }
>
>     return ans;
>   }
> }
>
> DEFUN_DLD (matern, args, ,
>    "Matern variogram function")
> {
>
>   /* Number of arguments*/
>   int nargs = args.length ();
>
>   /* Error checking*/
>   if (nargs!=3){
>     usage("matern (x, phi, kappa)");
>     return octave_value(0.0);
>   }
>
>   ColumnVector x = args(0).column_vector_value();
>   ColumnVector phi = args(1).column_vector_value();
>   ColumnVector kappa = args(2).column_vector_value();
>
>   int n = x.length();
>   ColumnVector retval = ColumnVector(n);
>   for( int i = 0; i < n; i++ ){
>     retval(i) = matern (x(i), phi(i), static_cast<int>(kappa(i)));
>   }
>
>   return retval;
> }
>
>
>
>



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------