outer product implementation
Muthiah Annamalai <[email protected]> Sun, 10 Sep 2006 12:59:51 -0500
| Newsgroups | gmane.comp.gnu.octave.sources |
|---|---|
| Message-ID | <1157911192.5512.20.camel@localhost> |
--=-FeeHII4f1E7Y0AidcIL8 Content-Type: text/plain Content-Transfer-Encoding: 7bit Paul, I have fixed the 2nd version of the outer product code as you had asked me to using C++ style inline functions etc & TexInfo documentation. The older code can be found at. http://www.cae.wisc.edu/pipermail/octave-sources/2006-April/000000.html Im attaching the new code & the same test case. Thanks Muthu --=-FeeHII4f1E7Y0AidcIL8 Content-Disposition: attachment; filename=outer.cpp Content-Type: text/x-c++src; name=outer.cpp; charset=UTF-8 Content-Transfer-Encoding: 7bit /* * (C) 2006, September, Muthiah Annamalai. <[email protected]> * An implementation of the 'outer'-product function as specified in the * octave-projects page, at www.octave.org. * * This program is free software; you can redistribute it and/or modify * it under the terms of the GNU General Public License as published by * the Free Software Foundation; either version 2 of the License, or * (at your option) any later version. * * This program is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with this program; if not, write to the Free Software * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA */ #include<iostream> #include<octave/oct.h> #include<octave/parse.h> #include<octave/dynamic-ld.h> #include<octave/oct-map.h> #include<octave/oct-stream.h> #include<octave/ov.h> #include<octave/parse.h> #include<octave/Matrix.h> #include<octave/ov-cx-mat.h> #include<octave/ov-list.h> #include<octave/ov-base-mat.h> #include<octave/quit.h> #include<octave/error.h> inline static bool is_vector(octave_value arg) { return ((arg.is_matrix_type()&& !arg.is_char_matrix() && (arg.columns()==1 || arg.rows()==1))? true : false); } DEFUN_DLD(outer,args,, "-*- texinfo -*-\n" "@deftypefn {Loadable Function} {@var{outer_product} =} outer (@var{x},@var{y})\n" "@deftypefnx {Loadable Function} {@var{outer_product} =} outer (@var{x},@var{y},@var{f})\n" "\n" "function outer(VecX,VecY,opt:function)\n" "Return a @var{M}x@var{N} matrix with @var{f}(@var{x},@var{y}) evaluated at x,y position\n" "for each @var{x},@var{y} in Vector @var{X},Vector @var{Y} respectively\n" "If function is not provided, use '*' as default operation\n" "Function can be provided using argument @var{f} \n" "using reference operator (@@) only.\n" "\n" "example: outer(@var{[1.0 2.0 3.0]},@var{[2.0 3.0 4.0]},@@gcd)\n" "compute the 3x3 matrix with each element in it being the\n" "GCD of the corresponding x and y elements of the matrix\n" "\n" "@end deftypefn\n" "@seealso{}") { int m=0,n=0,aLength=0; ColumnVector cv[2]; bool function_present=false; octave_function *func=NULL; octave_value_list erval; const char *func_name=NULL; bool is_commutative=false; aLength=args.length(); if(aLength < 2) { error("Usage: outer(VecX,VecY,opt:function)"); return erval; } if(!(is_vector(args(0)) || is_vector(args(1)))) { error("Usage: outer(VecX,VecY,opt:function);\n Invalid argument arg- [1 or 2]Expected Row/Column Vector"); return erval; } /* * FIXME: [Support for String-Name functions. Only references supported now] * Also allow functions to be passed by their names, and then replace the * name with a reference to the function-, later. * hope to use commutativity of the operators someday */ if(aLength>=3 && (args(2).is_function_handle() /* ||args(2).is_string() */)) { function_present=true; is_commutative=false; func=args(2).function_value(); } cv[0]=ColumnVector(args(0).vector_value()); cv[1]=ColumnVector(args(1).vector_value()); m=cv[0].length(); n=cv[1].length(); Matrix rmat(m,n); octave_value_list ovl; /* * If operator is commutative, then * do we want to handle this as a special * case? * */ if(function_present && func!=NULL) { /* try a dummy invoke on the function, and verify if the return value is there. Otherwise we report and error & die.*/ ovl(0)=0; ovl(1)=0; func_name = func->name().c_str(); if(feval(func,ovl).length() < 1) { error("The function %s must return atleast one value, on invocation\n",func_name); return octave_value(-1); } if(error_state) { error("Exception: Cannot find function specified"); return erval; } /* FIXME: <Vectorize Code> * Paul says do this thing as a single loop. * Iterate over Columns, taking one full column at a time, * and save the result on the columns of the matrix. * Atleast this is possible for builtins that take arguments * over columns/rows. */ for(int itr=0;itr<m;itr++) for(int icol=0;icol<n;icol++) { OCTAVE_QUIT; ovl(0)=cv[0].elem(itr); ovl(1)=cv[1].elem(icol); rmat.elem(itr,icol)=feval(func,ovl)(0).double_value(); if(error_state) { error("Exception: while evaluating %s with arguments %g %g\n",func_name,cv[0].elem(itr),cv[1].elem(icol)); return erval; } } } else { /* do a simple product */ rmat=cv[0]*cv[1].transpose(); } return octave_value(rmat); } /* Make command: g++ -fpic outer.cpp -shared -o outer.oct -Wall -ggdb -Wall -Wunused -Wconversion \ -fno-exceptions -DDEBUG=1 `mkoctfile -p INCFLAGS` */ --=-FeeHII4f1E7Y0AidcIL8 Content-Disposition: attachment; filename=outertest.m Content-Type: text/x-octave; name=outertest.m; charset=UTF-8 Content-Transfer-Encoding: 7bit #! /usr/bin/octave -q %outer disp('Help on Outer') help outer pause disp('Product- testing ;-) Imean this * product') outer([1 2],[1 2 3 4]) outer([1],[1 2 3]) outer([1 2 3],[1]) outer([1 2 3],[1 2 3]) %outer([1],[2]) %error this is caught. disp('Functions being Invoked') disp('Min') outer(rand(1,10),rand(1,10),@min) disp('IsEqual') outer([1 2],[1 2 3 4],@isequal) disp('Max') outer([1],[1 2 3],@max) disp('GCD ') outer([1.0 2.0 3.0],[2.0 3.0 4.0],@gcd) outer([14 20 35 40],[7],@gcd) disp('Square Matrix') outer([1:5],[1:5]) disp('One arg scalar other vector') outer(3,[1:5]) outer([1:5],3) disp('functions that dont return an argument! ') try outer([14 20 35 40],[7],@disp) %functions that dont return an argument! catch disp('Exception') end disp('non-existent functions throw errors at Octave prompt itself.'); try outer([14 20 35 40],[7],@bwand) %non-existent functions throw errors at Octave prompt itself. catch disp('Exception') end try disp('Both Scalar') outer([1],[2]) %error catch disp('Exception') end --=-FeeHII4f1E7Y0AidcIL8 Content-Type: text/plain; charset="us-ascii" MIME-Version: 1.0 Content-Transfer-Encoding: 7bit Content-Disposition: inline _______________________________________________ Octave-sources mailing list [email protected] https://www.cae.wisc.edu/mailman/listinfo/octave-sources --=-FeeHII4f1E7Y0AidcIL8--