|
6569
|
1 @node Dynamically Linked Functions |
|
|
2 @appendix Dynamically Linked Functions |
|
|
3 @cindex dynamic-linking |
|
|
4 |
|
|
5 Octave has the possibility of including compiled code as dynamically |
|
|
6 linked extensions and then using these extensions as if they were part |
|
|
7 of Octave itself. Octave has the option of directly calling C++ code |
|
|
8 through its native oct-file interface or C code through its mex |
|
|
9 interface. It can also indirectly call functions written in any other |
|
|
10 language through a simple wrapper. The reasons to write code in a |
|
|
11 compiled language might be either to link to an existing piece of code |
|
|
12 and allow it to be used within Octave, or to allow improved performance |
|
|
13 for key pieces of code. |
|
|
14 |
|
|
15 Before going further, you should first determine if you really need to |
|
|
16 use dynamically linked functions at all. Before proceeding with writing |
|
|
17 any dynamically linked function to improve performance you should |
|
|
18 address ask yourself |
|
|
19 |
|
|
20 @itemize @bullet |
|
|
21 @item |
|
|
22 Can I get the same functionality using the Octave scripting language only. |
|
|
23 @item |
|
|
24 Is it thoroughly optimized Octave code? Vectorization of Octave code, |
|
|
25 doesn't just make it concise, it generally significantly improves its |
|
|
26 performance. Above all, if loops must be used, make sure that the |
|
|
27 allocation of space for variables takes place outside the loops using an |
|
|
28 assignment to a like matrix or zeros. |
|
|
29 @item |
|
|
30 Does it make as much use as possible of existing built-in library |
|
|
31 routines? These are highly optimized and many do not carry the overhead |
|
|
32 of being interpreted. |
|
|
33 @item |
|
|
34 Does writing a dynamically linked function represent useful investment |
|
|
35 of your time, relative to staying in Octave? |
|
|
36 @end itemize |
|
|
37 |
|
|
38 Also, as oct- and mex-files are dynamically linked to octave, they |
|
|
39 introduce to possibility of having Octave abort due to coding errors in |
|
|
40 the user code. For example a segmentation violation in the users code |
|
|
41 will cause Octave to abort. |
|
|
42 |
|
|
43 @menu |
|
|
44 * Oct-Files:: |
|
|
45 * Mex-Files:: |
|
|
46 * Standalone Programs:: |
|
|
47 @end menu |
|
|
48 |
|
|
49 @node Oct-Files |
|
|
50 @section Oct-Files |
|
|
51 @cindex oct-files |
|
|
52 @cindex mkoctfile |
|
|
53 @cindex oct |
|
|
54 |
|
|
55 @menu |
|
|
56 * Getting Started with Oct-Files:: |
|
|
57 * Matrices and Arrays in Oct-Files:: |
|
|
58 * Using Sparse Matrices in Oct-Files:: |
|
|
59 * Using Strings in Oct-Files:: |
|
|
60 * Cell Arrays in Oct-Files:: |
|
|
61 * Using Structures in Oct-Files:: |
|
|
62 * Accessing Global Variables in Oct-Files:: |
|
|
63 * Calling Octave Functions from Oct-Files:: |
|
|
64 * Calling External Code from Oct-Files:: |
|
|
65 * Allocating Local Memory in Oct-Files:: |
|
|
66 * Input Parameter Checking in Oct-Files:: |
|
|
67 * Exception and Error Handling in Oct-Files:: |
|
|
68 * Documentation and Test of Oct-Files:: |
|
|
69 * Application Programming Interface for Oct-Files:: |
|
|
70 @end menu |
|
|
71 |
|
|
72 @node Getting Started with Oct-Files |
|
|
73 @subsection Getting Started with Oct-Files |
|
|
74 |
|
|
75 The basic command to build oct-files is @code{mkoctfile} and it can be |
|
|
76 call from within octave or from the command line. |
|
|
77 |
|
|
78 @DOCSTRING(mkoctfile) |
|
|
79 |
|
|
80 Consider the short example |
|
|
81 |
|
|
82 @example |
|
|
83 @group |
|
|
84 #include <octave/oct.h> |
|
|
85 |
|
|
86 DEFUN_DLD (helloworld, args, nargout, |
|
|
87 "Hello World Help String") |
|
|
88 @{ |
|
|
89 int nargin = args.length (); |
|
|
90 octave_stdout << "Hello World has " << nargin |
|
|
91 << " input arguments and " |
|
|
92 << nargout << " output arguments.\n"; |
|
|
93 return octave_value_list (); |
|
|
94 @} |
|
|
95 @end group |
|
|
96 @end example |
|
|
97 |
|
|
98 This example although short introduces the basics of writing a C++ |
|
|
99 function that can be dynamically linked to Octave. The easiest way to |
|
|
100 make available most of the definitions that might be necessary for an |
|
|
101 oct-file in Octave is to use the @code{#include <octave/oct.h>} |
|
|
102 header. |
|
|
103 |
|
|
104 The macro that defines the entry point into the dynamically loaded |
|
|
105 function is DEFUN_DLD. This macro takes four arguments, these being |
|
|
106 |
|
|
107 @enumerate 1 |
|
|
108 @item The function name as it will be seen in Octave, |
|
|
109 @item The list of arguments to the function of type octave_value_list, |
|
|
110 @item The number of output arguments, which can and often is omitted if |
|
|
111 not used, and |
|
|
112 @item The string that will be seen as the help text of the function. |
|
|
113 @end enumerate |
|
|
114 |
|
|
115 The return type of functions defined with DEFUN_DLD is always |
|
|
116 octave_value_list. |
|
|
117 |
|
|
118 There are a couple of important considerations in the choice of function |
|
|
119 name. Firstly, it must be a valid Octave function name and so must be a |
|
|
120 sequence of letters, digits and underscores, not starting with a |
|
|
121 digit. Secondly, as Octave uses the function name to define the filename |
|
|
122 it attempts to find the function in, the function name in the DEFUN_DLD |
|
|
123 macro must match the filename of the oct-file. Therefore, the above |
|
|
124 function should be in a file helloworld.cc, and it should be compiled to |
|
|
125 an oct-file using the command |
|
|
126 |
|
|
127 @example |
|
|
128 mkoctfile helloworld.cc |
|
|
129 @end example |
|
|
130 |
|
|
131 This will create a file call helloworld.oct, that is the compiled |
|
|
132 version of the function. It should be noted that it is perfectly |
|
|
133 acceptable to have more than one DEFUN_DLD function in a source |
|
|
134 file. However, there must either be a symbolic link to the oct-file for |
|
|
135 each of the functions defined in the source code with the DEFUN_DLD |
|
|
136 macro or the autoload (@ref{Function Files}) function should be used. |
|
|
137 |
|
|
138 The rest of this function then shows how to find the number of input |
|
|
139 arguments, how to print through the octave pager, and return from the |
|
|
140 function. After compiling this function as above, an example of its use |
|
|
141 is |
|
|
142 |
|
|
143 @example |
|
|
144 @group |
|
|
145 helloworld(1,2,3) |
|
|
146 @result{} Hello World has 3 input arguments and 0 output arguments. |
|
|
147 @end group |
|
|
148 @end example |
|
|
149 |
|
|
150 @node Matrices and Arrays in Oct-Files |
|
|
151 @subsection Matrices and Arrays in Oct-Files |
|
|
152 |
|
|
153 Octave supports a number of different array and matrix classes, the |
|
|
154 majority of which are based on the Array class. The exception is the |
|
|
155 sparse matrix types discussed separately below. There are three basic |
|
|
156 matrix types |
|
|
157 |
|
|
158 @table @asis |
|
|
159 @item Matrix |
|
|
160 A double precision matrix class defined in dMatrix.h, |
|
|
161 @item ComplexMatrix |
|
|
162 A complex matrix class defined in CMatrix.h, and |
|
|
163 @item BoolMatrix |
|
|
164 A boolean matrix class defined in boolMatrix.h. |
|
|
165 @end table |
|
|
166 |
|
|
167 These are the basic two-dimensional matrix types of octave. In |
|
|
168 additional there are a number of multi-dimensional array types, these |
|
|
169 being |
|
|
170 |
|
|
171 @table @asis |
|
|
172 @item NDArray |
|
|
173 A double precision array class defined in dNDArray.h, |
|
|
174 @item ComplexNDarray |
|
|
175 A complex array class defined in CNDArray.h, |
|
|
176 @item boolNDArray |
|
|
177 A boolean array class defined in boolNDArray.h, |
|
|
178 @item int8NDArray, int16NDArray, int32NDArray, int64NDArray |
|
|
179 8, 16, 32 and 64-bit signed array classes defined in int8NDArray.h, etc, and |
|
|
180 @item uint8NDArray, uint16NDArray, uint32NDArray, uint64NDArray |
|
|
181 8, 16, 32 and 64-bit unsigned array classes defined in uint8NDArray.h, |
|
|
182 etc. |
|
|
183 @end table |
|
|
184 |
|
|
185 There are several basic means of constructing matrices of |
|
|
186 multi-dimensional arrays. Considering the Matrix type as an example |
|
|
187 |
|
|
188 @itemize @bullet |
|
|
189 @item |
|
|
190 We can create an empty matrix or array with the empty constructor. For |
|
|
191 example |
|
|
192 |
|
|
193 @example |
|
|
194 Matrix a; |
|
|
195 @end example |
|
|
196 |
|
|
197 This can be used on all matrix and array types |
|
|
198 @item |
|
|
199 Define the dimensions of the matrix or array with a dim_vector. For |
|
|
200 example |
|
|
201 |
|
|
202 @example |
|
|
203 @group |
|
|
204 dim_vector dv(2); |
|
|
205 dv(0) = 2; dv(1) = 2; |
|
|
206 Matrix a(dv); |
|
|
207 @end group |
|
|
208 @end example |
|
|
209 |
|
|
210 This can be used on all matrix and array types |
|
|
211 @item |
|
|
212 Define the number of rows and columns in the Matrix. For example |
|
|
213 |
|
|
214 @example |
|
|
215 Matrix a(2,2) |
|
|
216 @end example |
|
|
217 |
|
|
218 However, this constructor can only be used with the matrix types. |
|
|
219 @end itemize |
|
|
220 |
|
|
221 These types all share a number of basic methods and operators, a |
|
|
222 selection of which include |
|
|
223 |
|
|
224 |
|
|
225 @table @asis |
|
|
226 @item T& operator (octave_idx_type), T& elem(octave_idx_type) |
|
|
227 The () operator or elem method allow the values of the matrix or array |
|
|
228 to be read or set. These can take a single argument, which is of type |
|
|
229 octave_idx_type, that is the index into the matrix or |
|
|
230 array. Additionally, the matrix type allows two argument versions of the |
|
|
231 () operator and elem method, giving the row and column index of the |
|
|
232 value to obtain or set. |
|
|
233 |
|
|
234 Note that these function do significant error checking and so in some |
|
|
235 circumstances the user might prefer the access the data of the array or |
|
|
236 matrix directly through the fortran_vec method discussed below. |
|
|
237 @item octave_idx_type nelem () |
|
|
238 The total number of elements in the matrix or array. |
|
|
239 @item size_t byte_size () |
|
|
240 The number of bytes used to store the matrix or array. |
|
|
241 @item dim_vector dims() |
|
|
242 The dimensions of the matrix or array in value of type dim_vector. |
|
|
243 @item void resize (const dim_vector&) |
|
|
244 A method taking either an argument of type dim_vector, or in the case of |
|
|
245 a matrix two arguments of type octave_idx_type defining the number of |
|
|
246 rows and columns in the matrix. |
|
|
247 @item T* fortran_vec () |
|
|
248 This method returns a pointer to the underlying data of the matrix or a |
|
|
249 array so that it can be manipulated directly, either within Octave or by |
|
|
250 an external library. |
|
|
251 @end table |
|
|
252 |
|
|
253 Operators such an +, -, or * can be used on the majority of the above |
|
|
254 types. In addition there are a number of methods that are of interest |
|
|
255 only for matrices such as transpose, hermitian, solve, etc. |
|
|
256 |
|
|
257 The typical way to extract a matrix or array from the input arguments of |
|
|
258 DEFUN_DLD function is as follows |
|
|
259 |
|
|
260 @example |
|
|
261 @group |
|
|
262 #include <octave/oct.h> |
|
|
263 |
|
|
264 DEFUN_DLD (addtwomatrices, args, , "Add A to B") |
|
|
265 @{ |
|
|
266 int nargin = args.length (); |
|
|
267 if (nargin != 2) |
|
|
268 print_usage (); |
|
|
269 else |
|
|
270 @{ |
|
|
271 NDArray A = args(0).array_value(); |
|
|
272 NDArray B = args(1).array_value(); |
|
|
273 if (! error_state) |
|
|
274 return octave_value(A + B); |
|
|
275 @} |
|
|
276 return octave_value_list (); |
|
|
277 @} |
|
|
278 @end group |
|
|
279 @end example |
|
|
280 |
|
|
281 To avoid segmentation faults causing Octave to abort, this function |
|
|
282 explicitly checks that there are sufficient arguments available before |
|
|
283 accessing these arguments. It then obtains two multi-dimensional arrays |
|
|
284 of type NDArray and adds these together. Note that the array_value |
|
|
285 method is called without using the is_matrix_type type, and instead the |
|
|
286 error_state is checked before returning @code{A + B}. The reason to |
|
|
287 prefer this is that the arguments might be a type that is not an |
|
|
288 NDArray, but it would make sense to convert it to one. The array_value |
|
|
289 method allows this conversion to be performed transparently if possible, |
|
|
290 and sets error_state if it is not. |
|
|
291 |
|
|
292 @code{A + B}, operating on two NDArray's returns an NDArray, which is |
|
|
293 cast to an octave_value on the return from the function. An example of |
|
|
294 the use of this demonstration function is |
|
|
295 |
|
|
296 @example |
|
|
297 @group |
|
|
298 addtwomatrices(ones(2,2),ones(2,2)) |
|
|
299 @result{} 2 2 |
|
|
300 2 2 |
|
|
301 @end group |
|
|
302 @end example |
|
|
303 |
|
|
304 A list of the basic Matrix and Array types, the methods to extract these |
|
|
305 from an octave_value and the associated header is listed below. |
|
|
306 |
|
|
307 @multitable @columnfractions .3 .4 .3 |
|
|
308 @item RowVector @tab row_vector_value @tab dRowVector.h |
|
|
309 @item ComplexRowVector @tab complex_row_vector_value @tab CRowVector.h |
|
|
310 @item ColumnVector @tab column_vector_value @tab dColVector.h |
|
|
311 @item ComplexColumnVector @tab complex_column_vector_value @tab CColVector.h |
|
|
312 @item Matrix @tab matrix_value @tab dMatrix.h |
|
|
313 @item ComplexMatrix @tab complex_matrix_value @tab CMatrix.h |
|
|
314 @item boolMatrix @tab bool_matrix_value @tab boolMatrix.h |
|
|
315 @item charMatrix @tab char_matrix_value @tab chMatrix.h |
|
|
316 @item NDArray @tab array_value @tab dNDArray.h |
|
|
317 @item ComplexNDArray @tab complex_array_value @tab CNDArray.h |
|
|
318 @item boolNDArray @tab bool_array_value @tab boolNDArray.h |
|
|
319 @item charNDArray @tab char_array_value @tab charNDArray.h |
|
|
320 @item int8NDArray @tab int8_array_value @tab int8NDArray.h |
|
|
321 @item int16NDArray @tab int16_array_value @tab int16NDArray.h |
|
|
322 @item int32NDArray @tab int32_array_value @tab int32NDArray.h |
|
|
323 @item int64NDArray @tab int64_array_value @tab int64NDArray.h |
|
|
324 @item uint8NDArray @tab uint8_array_value @tab uint8NDArray.h |
|
|
325 @item uint16NDArray @tab uint16_array_value @tab uint16NDArray.h |
|
|
326 @item uint32NDArray @tab uint32_array_value @tab uint32NDArray.h |
|
|
327 @item uint64NDArray @tab uint64_array_value @tab uint64NDArray.h |
|
|
328 @end multitable |
|
|
329 |
|
|
330 @node Using Sparse Matrices in Oct-Files |
|
|
331 @subsection Using Sparse Matrices in Oct-Files |
|
|
332 |
|
|
333 There are three classes of sparse objects that are of interest to the |
|
|
334 user. |
|
|
335 |
|
|
336 @table @asis |
|
|
337 @item SparseMatrix |
|
|
338 A double precision sparse matrix class |
|
|
339 @item SparseComplexMatrix |
|
|
340 A complex sparse matrix class |
|
|
341 @item SparseBoolMatrix |
|
|
342 A boolean sparse matrix class |
|
|
343 @end table |
|
|
344 |
|
|
345 All of these classes inherit from the @code{Sparse<T>} template class, |
|
|
346 and so all have similar capabilities and usage. The @code{Sparse<T>} |
|
|
347 class was based on Octave @code{Array<T>} class, and so users familiar |
|
|
348 with Octave's Array classes will be comfortable with the use of |
|
|
349 the sparse classes. |
|
|
350 |
|
|
351 The sparse classes will not be entirely described in this section, due |
|
|
352 to their similar with the existing Array classes. However, there are a |
|
|
353 few differences due the different nature of sparse objects, and these |
|
|
354 will be described. Firstly, although it is fundamentally possible to |
|
|
355 have N-dimensional sparse objects, the Octave sparse classes do |
|
|
356 not allow them at this time. So all operations of the sparse classes |
|
|
357 must be 2-dimensional. This means that in fact @code{SparseMatrix} is |
|
|
358 similar to Octave's @code{Matrix} class rather than its |
|
|
359 @code{NDArray} class. |
|
|
360 |
|
|
361 @menu |
|
|
362 * OctDifferences:: The Differences between the Array and Sparse Classes |
|
|
363 * OctCreation:: Creating Spare Matrices in Oct-Files |
|
|
364 * OctUse:: Using Sparse Matrices in Oct-Files |
|
|
365 @end menu |
|
|
366 |
|
|
367 @node OctDifferences |
|
|
368 @subsubsection The Differences between the Array and Sparse Classes |
|
|
369 |
|
|
370 The number of elements in a sparse matrix is considered to be the number |
|
|
371 of non-zero elements rather than the product of the dimensions. Therefore |
|
|
372 |
|
|
373 @example |
|
|
374 SparseMatrix sm; |
|
|
375 @dots{} |
|
|
376 int nel = sm.nelem (); |
|
|
377 @end example |
|
|
378 |
|
|
379 returns the number of non-zero elements. If the user really requires the |
|
|
380 number of elements in the matrix, including the non-zero elements, they |
|
|
381 should use @code{numel} rather than @code{nelem}. Note that for very |
|
|
382 large matrices, where the product of the two dimensions is large that |
|
|
383 the representation of the an unsigned int, then @code{numel} can overflow. |
|
|
384 An example is @code{speye(1e6)} which will create a matrix with a million |
|
|
385 rows and columns, but only a million non-zero elements. Therefore the |
|
|
386 number of rows by the number of columns in this case is more than two |
|
|
387 hundred times the maximum value that can be represented by an unsigned int. |
|
|
388 The use of @code{numel} should therefore be avoided useless it is known |
|
|
389 it won't overflow. |
|
|
390 |
|
|
391 Extreme care must be take with the elem method and the "()" operator, |
|
|
392 which perform basically the same function. The reason is that if a |
|
|
393 sparse object is non-const, then Octave will assume that a |
|
|
394 request for a zero element in a sparse matrix is in fact a request |
|
|
395 to create this element so it can be filled. Therefore a piece of |
|
|
396 code like |
|
|
397 |
|
|
398 @example |
|
|
399 SparseMatrix sm; |
|
|
400 @dots{} |
|
|
401 for (int j = 0; j < nc; j++) |
|
|
402 for (int i = 0; i < nr; i++) |
|
|
403 std::cerr << " (" << i << "," << j << "): " << sm(i,j) |
|
|
404 << std::endl; |
|
|
405 @end example |
|
|
406 |
|
|
407 is a great way of turning the sparse matrix into a dense one, and a |
|
|
408 very slow way at that since it reallocates the sparse object at each |
|
|
409 zero element in the matrix. |
|
|
410 |
|
|
411 An easy way of preventing the above from happening is to create a temporary |
|
|
412 constant version of the sparse matrix. Note that only the container for |
|
|
413 the sparse matrix will be copied, while the actual representation of the |
|
|
414 data will be shared between the two versions of the sparse matrix. So this |
|
|
415 is not a costly operation. For example, the above would become |
|
|
416 |
|
|
417 @example |
|
|
418 SparseMatrix sm; |
|
|
419 @dots{} |
|
|
420 const SparseMatrix tmp (sm); |
|
|
421 for (int j = 0; j < nc; j++) |
|
|
422 for (int i = 0; i < nr; i++) |
|
|
423 std::cerr << " (" << i << "," << j << "): " << tmp(i,j) |
|
|
424 << std::endl; |
|
|
425 @end example |
|
|
426 |
|
|
427 Finally, as the sparse types aren't just represented as a contiguous |
|
|
428 block of memory, the @code{fortran_vec} method of the @code{Array<T>} |
|
|
429 is not available. It is however replaced by three separate methods |
|
|
430 @code{ridx}, @code{cidx} and @code{data}, that access the raw compressed |
|
|
431 column format that the Octave sparse matrices are stored in. |
|
|
432 Additionally, these methods can be used in a manner similar to @code{elem}, |
|
|
433 to allow the matrix to be accessed or filled. However, in that case it is |
|
|
434 up to the user to respect the sparse matrix compressed column format |
|
|
435 discussed previous. |
|
|
436 |
|
|
437 @node OctCreation |
|
|
438 @subsubsection Creating Spare Matrices in Oct-Files |
|
|
439 |
|
|
440 The user has several alternatives in how to create a sparse matrix. |
|
|
441 They can first create the data as three vectors representing the |
|
|
442 row and column indexes and the data, and from those create the matrix. |
|
|
443 Or alternatively, they can create a sparse matrix with the appropriate |
|
|
444 amount of space and then fill in the values. Both techniques have their |
|
|
445 advantages and disadvantages. |
|
|
446 |
|
|
447 An example of how to create a small sparse matrix with the first technique |
|
|
448 might be seen the example |
|
|
449 |
|
|
450 @example |
|
|
451 int nz = 4, nr = 3, nc = 4; |
|
|
452 ColumnVector ridx (nz); |
|
|
453 ColumnVector cidx (nz); |
|
|
454 ColumnVector data (nz); |
|
|
455 |
|
|
456 ridx(0) = 0; ridx(1) = 0; ridx(2) = 1; ridx(3) = 2; |
|
|
457 cidx(0) = 0; cidx(1) = 1; cidx(2) = 3; cidx(3) = 3; |
|
|
458 data(0) = 1; data(1) = 2; data(2) = 3; data(3) = 4; |
|
|
459 |
|
|
460 SparseMatrix sm (data, ridx, cidx, nr, nc); |
|
|
461 @end example |
|
|
462 |
|
|
463 which creates the matrix given in section @ref{Storage}. Note that |
|
|
464 the compressed matrix format is not used at the time of the creation |
|
|
465 of the matrix itself, however it is used internally. |
|
|
466 |
|
|
467 As previously mentioned, the values of the sparse matrix are stored |
|
|
468 in increasing column-major ordering. Although the data passed by the |
|
|
469 user does not need to respect this requirement, the pre-sorting the |
|
|
470 data significantly speeds up the creation of the sparse matrix. |
|
|
471 |
|
|
472 The disadvantage of this technique of creating a sparse matrix is |
|
|
473 that there is a brief time where two copies of the data exists. Therefore |
|
|
474 for extremely memory constrained problems this might not be the right |
|
|
475 technique to create the sparse matrix. |
|
|
476 |
|
|
477 The alternative is to first create the sparse matrix with the desired |
|
|
478 number of non-zero elements and then later fill those elements in. The |
|
|
479 easiest way to do this is |
|
|
480 |
|
|
481 @example |
|
|
482 int nz = 4, nr = 3, nc = 4; |
|
|
483 SparseMatrix sm (nr, nc, nz); |
|
|
484 sm(0,0) = 1; sm(0,1) = 2; sm(1,3) = 3; sm(2,3) = 4; |
|
|
485 @end example |
|
|
486 |
|
|
487 That creates the same matrix as previously. Again, although it is not |
|
|
488 strictly necessary, it is significantly faster if the sparse matrix is |
|
|
489 created in this manner that the elements are added in column-major |
|
|
490 ordering. The reason for this is that if the elements are inserted |
|
|
491 at the end of the current list of known elements then no element |
|
|
492 in the matrix needs to be moved to allow the new element to be |
|
|
493 inserted. Only the column indexes need to be updated. |
|
|
494 |
|
|
495 There are a few further points to note about this technique of creating |
|
|
496 a sparse matrix. Firstly, it is not illegal to create a sparse matrix |
|
|
497 with fewer elements than are actually inserted in the matrix. Therefore |
|
|
498 |
|
|
499 @example |
|
|
500 int nz = 4, nr = 3, nc = 4; |
|
|
501 SparseMatrix sm (nr, nc, 0); |
|
|
502 sm(0,0) = 1; sm(0,1) = 2; sm(1,3) = 3; sm(2,3) = 4; |
|
|
503 @end example |
|
|
504 |
|
|
505 is perfectly legal. However it is a very bad idea. The reason is that |
|
|
506 as each new element is added to the sparse matrix the space allocated |
|
|
507 to it is increased by reallocating the memory. This is an expensive |
|
|
508 operation, that will significantly slow this means of creating a sparse |
|
|
509 matrix. Furthermore, it is not illegal to create a sparse matrix with |
|
|
510 too much storage, so having @var{nz} above equaling 6 is also legal. |
|
|
511 The disadvantage is that the matrix occupies more memory than strictly |
|
|
512 needed. |
|
|
513 |
|
|
514 It is not always easy to known the number of non-zero elements prior |
|
|
515 to filling a matrix. For this reason the additional storage for the |
|
|
516 sparse matrix can be removed after its creation with the |
|
|
517 @dfn{maybe_compress} function. Furthermore, the maybe_compress can |
|
|
518 deallocate the unused storage, but it can equally remove zero elements |
|
|
519 from the matrix. The removal of zero elements from the matrix is |
|
|
520 controlled by setting the argument of the @dfn{maybe_compress} function |
|
|
521 to be 'true'. However, the cost of removing the zeros is high because it |
|
|
522 implies resorting the elements. Therefore, if possible it is better |
|
|
523 is the user doesn't add the zeros in the first place. An example of |
|
|
524 the use of @dfn{maybe_compress} is |
|
|
525 |
|
|
526 @example |
|
|
527 int nz = 6, nr = 3, nc = 4; |
|
|
528 SparseMatrix sm1 (nr, nc, nz); |
|
|
529 sm1(0,0) = 1; sm1(0,1) = 2; sm1(1,3) = 3; sm1(2,3) = 4; |
|
|
530 sm1.maybe_compress (); // No zero elements were added |
|
|
531 |
|
|
532 SparseMatrix sm2 (nr, nc, nz); |
|
|
533 sm2(0,0) = 1; sm2(0,1) = 2; sm(0,2) = 0; sm(1,2) = 0; |
|
|
534 sm1(1,3) = 3; sm1(2,3) = 4; |
|
|
535 sm2.maybe_compress (true); // Zero elements were added |
|
|
536 @end example |
|
|
537 |
|
|
538 The use of the @dfn{maybe_compress} function should be avoided if |
|
|
539 possible, as it will slow the creation of the matrices. |
|
|
540 |
|
|
541 A third means of creating a sparse matrix is to work directly with |
|
|
542 the data in compressed row format. An example of this technique might |
|
|
543 be |
|
|
544 |
|
|
545 @c Note the @verbatim environment is a relatively new addition to texinfo. |
|
|
546 @c Therefore use the @example environment and replace @, with @@, |
|
|
547 @c { with @{, etc |
|
|
548 |
|
|
549 @example |
|
|
550 octave_value arg; |
|
|
551 |
|
|
552 @dots{} |
|
|
553 |
|
|
554 int nz = 6, nr = 3, nc = 4; // Assume we know the max no nz |
|
|
555 SparseMatrix sm (nr, nc, nz); |
|
|
556 Matrix m = arg.matrix_value (); |
|
|
557 |
|
|
558 int ii = 0; |
|
|
559 sm.cidx (0) = 0; |
|
|
560 for (int j = 1; j < nc; j++) |
|
|
561 @{ |
|
|
562 for (int i = 0; i < nr; i++) |
|
|
563 @{ |
|
|
564 double tmp = foo (m(i,j)); |
|
|
565 if (tmp != 0.) |
|
|
566 @{ |
|
|
567 sm.data(ii) = tmp; |
|
|
568 sm.ridx(ii) = i; |
|
|
569 ii++; |
|
|
570 @} |
|
|
571 @} |
|
|
572 sm.cidx(j+1) = ii; |
|
|
573 @} |
|
|
574 sm.maybe_compress (); // If don't know a-priori the final no of nz. |
|
|
575 @end example |
|
|
576 |
|
|
577 which is probably the most efficient means of creating the sparse matrix. |
|
|
578 |
|
|
579 Finally, it might sometimes arise that the amount of storage initially |
|
|
580 created is insufficient to completely store the sparse matrix. Therefore, |
|
|
581 the method @code{change_capacity} exists to reallocate the sparse memory. |
|
|
582 The above example would then be modified as |
|
|
583 |
|
|
584 @example |
|
|
585 octave_value arg; |
|
|
586 |
|
|
587 @dots{} |
|
|
588 |
|
|
589 int nz = 6, nr = 3, nc = 4; // Assume we know the max no nz |
|
|
590 SparseMatrix sm (nr, nc, nz); |
|
|
591 Matrix m = arg.matrix_value (); |
|
|
592 |
|
|
593 int ii = 0; |
|
|
594 sm.cidx (0) = 0; |
|
|
595 for (int j = 1; j < nc; j++) |
|
|
596 @{ |
|
|
597 for (int i = 0; i < nr; i++) |
|
|
598 @{ |
|
|
599 double tmp = foo (m(i,j)); |
|
|
600 if (tmp != 0.) |
|
|
601 @{ |
|
|
602 if (ii == nz) |
|
|
603 @{ |
|
|
604 nz += 2; // Add 2 more elements |
|
|
605 sm.change_capacity (nz); |
|
|
606 @} |
|
|
607 sm.data(ii) = tmp; |
|
|
608 sm.ridx(ii) = i; |
|
|
609 ii++; |
|
|
610 @} |
|
|
611 @} |
|
|
612 sm.cidx(j+1) = ii; |
|
|
613 @} |
|
|
614 sm.maybe_mutate (); // If don't know a-priori the final no of nz. |
|
|
615 @end example |
|
|
616 |
|
|
617 Note that both increasing and decreasing the number of non-zero elements in |
|
|
618 a sparse matrix is expensive, as it involves memory reallocation. Also as |
|
|
619 parts of the matrix, though not its entirety, exist as the old and new copy |
|
|
620 at the same time, additional memory is needed. Therefore if possible this |
|
|
621 should be avoided. |
|
|
622 |
|
|
623 @node OctUse |
|
|
624 @subsubsection Using Sparse Matrices in Oct-Files |
|
|
625 |
|
|
626 Most of the same operators and functions on sparse matrices that are |
|
|
627 available from the Octave are equally available with oct-files. |
|
|
628 The basic means of extracting a sparse matrix from an @code{octave_value} |
|
|
629 and returning them as an @code{octave_value}, can be seen in the |
|
|
630 following example |
|
|
631 |
|
|
632 @example |
|
|
633 octave_value_list retval; |
|
|
634 |
|
|
635 SparseMatrix sm = args(0).sparse_matrix_value (); |
|
|
636 SparseComplexMatrix scm = args(1).sparse_complex_matrix_value (); |
|
|
637 SparseBoolMatrix sbm = args(2).sparse_bool_matrix_value (); |
|
|
638 |
|
|
639 @dots{} |
|
|
640 |
|
|
641 retval(2) = sbm; |
|
|
642 retval(1) = scm; |
|
|
643 retval(0) = sm; |
|
|
644 @end example |
|
|
645 |
|
|
646 The conversion to an octave-value is handled by the sparse |
|
|
647 @code{octave_value} constructors, and so no special care is needed. |
|
|
648 |
|
|
649 @node Using Strings in Oct-Files |
|
|
650 @subsection Using Strings in Oct-Files |
|
|
651 |
|
|
652 In Octave a string is just a special Array class. Consider the example |
|
|
653 |
|
|
654 @example |
|
|
655 @group |
|
|
656 #include <octave/oct.h> |
|
|
657 |
|
|
658 DEFUN_DLD (stringdemo, args, , "String Demo") |
|
|
659 @{ |
|
|
660 int nargin = args.length(); |
|
|
661 octave_value_list retval; |
|
|
662 |
|
|
663 if (nargin != 1) |
|
|
664 print_usage (); |
|
|
665 else |
|
|
666 @{ |
|
|
667 charMatrix ch = args(0).char_matrix_value (); |
|
|
668 |
|
|
669 if (! error_state) |
|
|
670 @{ |
|
|
671 if (args(0).is_sq_string ()) |
|
|
672 retval(1) = octave_value (ch, true); |
|
|
673 else |
|
|
674 retval(1) = octave_value (ch, true, '\''); |
|
|
675 |
|
|
676 octave_idx_type nr = ch.rows(); |
|
|
677 for (octave_idx_type i = 0; i < nr / 2; i++) |
|
|
678 @{ |
|
|
679 std::string tmp = ch.row_as_string (i); |
|
|
680 ch.insert (ch.row_as_string(nr - i - 1).c_str(), i, 0); |
|
|
681 ch.insert (tmp.c_str(), nr - i - 1, 0); |
|
|
682 @} |
|
|
683 retval(0) = octave_value (ch, true); |
|
|
684 @} |
|
|
685 @} |
|
|
686 return retval; |
|
|
687 @} |
|
|
688 @end group |
|
|
689 @end example |
|
|
690 |
|
|
691 An example of the of the use of this function is |
|
|
692 |
|
|
693 @example |
|
|
694 @group |
|
|
695 s0 = ["First String";"Second String"]; |
|
|
696 [s1,s2] = stringdemo (s0) |
|
|
697 @result{}s1 = Second String |
|
|
698 First String |
|
|
699 |
|
|
700 @result{}s2 = First String |
|
|
701 Second String |
|
|
702 |
|
|
703 typeinfo (s2) |
|
|
704 @result{} sq_string |
|
|
705 typeinfo (s1) |
|
|
706 @result{} string |
|
|
707 @end group |
|
|
708 @end example |
|
|
709 |
|
|
710 One additional complication of strings in Octave is the difference |
|
|
711 between single quoted and double quoted strings. To find out if an |
|
|
712 octave_value contains a single or double quoted string an example is |
|
|
713 |
|
|
714 @example |
|
|
715 @group |
|
|
716 if (args(0).is_sq_string()) |
|
|
717 octave_stdout << "First argument is a singularly quoted string\n"; |
|
|
718 else if (args(0).is_dq_string()) |
|
|
719 octave_stdout << "First argument is a doubly quoted string\n"; |
|
|
720 @end group |
|
|
721 @end example |
|
|
722 |
|
|
723 Note however, that both types of strings are represented by the |
|
|
724 charNDArray type, and so when assigning to an octave_value, the type of |
|
|
725 string should be specified. For example |
|
|
726 |
|
|
727 @example |
|
|
728 @group |
|
|
729 octave_value_list retval; |
|
|
730 charNDArray c; |
|
|
731 @dots{} |
|
|
732 retval(0) = octave_value (ch, true); // Create a double quoted string |
|
|
733 retval(1) = octave_value (ch, true, "'"); // Create single quoted string |
|
|
734 @end group |
|
|
735 @end example |
|
|
736 |
|
|
737 @node Cell Arrays in Oct-Files |
|
|
738 @subsection Cell Arrays in Oct-Files |
|
|
739 |
|
|
740 Octave's cell type is equally accessible within an oct-files. A cell |
|
|
741 array is just an array of octave_values, and so each element of the cell |
|
|
742 array can then be treated just like any other octave_value. A simple |
|
|
743 example is |
|
|
744 |
|
|
745 @example |
|
|
746 @group |
|
|
747 #include <octave/oct.h> |
|
|
748 #include <octave/Cell.h> |
|
|
749 |
|
|
750 DEFUN_DLD (celldemo, args, , "Cell Demo") |
|
|
751 @{ |
|
|
752 octave_value_list retval; |
|
|
753 int nargin = args.length(); |
|
|
754 |
|
|
755 if (nargin != 1) |
|
|
756 print_usage (); |
|
|
757 else |
|
|
758 @{ |
|
|
759 Cell c = args (0).cell_value (); |
|
|
760 if (! error_state) |
|
|
761 for (octave_idx_type i = 0; i < c.nelem (); i++) |
|
|
762 retval(i) = c.elem (i); |
|
|
763 @} |
|
|
764 return retval; |
|
|
765 @} |
|
|
766 @end group |
|
|
767 @end example |
|
|
768 |
|
|
769 Note that cell arrays are used less often in standard oct-files and so |
|
|
770 the Cell.h header file must be explicitly included. The rest of this |
|
|
771 example extracts the octave_values one by one from the cell array and |
|
|
772 returns be as individual return arguments. For example consider |
|
|
773 |
|
|
774 @example |
|
|
775 @group |
|
|
776 [b1,b2,b3] = celldemo(@{1, [1,2], "test"@}) |
|
|
777 @result{} |
|
|
778 b1 = 1 |
|
|
779 b2 = |
|
|
780 |
|
|
781 1 2 |
|
|
782 |
|
|
783 b3 = test |
|
|
784 @end group |
|
|
785 @end example |
|
|
786 |
|
|
787 @node Using Structures in Oct-Files |
|
|
788 @subsection Using Structures in Oct-Files |
|
|
789 |
|
|
790 A structure in Octave is map between a number of fields represented and |
|
|
791 their values. The "Standard Template Library" map class is used, with |
|
|
792 the pair consisting of a std::string and an octave Cell variable. |
|
|
793 |
|
|
794 A simple example demonstrating the use of structures within oct-files is |
|
|
795 |
|
|
796 @example |
|
|
797 @group |
|
|
798 #include <octave/oct.h> |
|
|
799 #include <octave/ov-struct.h> |
|
|
800 |
|
|
801 DEFUN_DLD (structdemo, args, , "Struct demo.") |
|
|
802 @{ |
|
|
803 int nargin = args.length (); |
|
|
804 octave_value retval; |
|
|
805 |
|
|
806 if (nargin != 2) |
|
|
807 print_usage (); |
|
|
808 else |
|
|
809 @{ |
|
|
810 Octave_map arg0 = args(0).map_value (); |
|
|
811 std::string arg1 = args(1).string_value (); |
|
|
812 |
|
|
813 if (! error_state && arg0.contains (arg1)) |
|
|
814 @{ |
|
|
815 // The following two lines might be written as |
|
|
816 // octave_value tmp; |
|
|
817 // for (Octave_map::iterator p0 = arg0.begin() ; |
|
|
818 // p0 != arg0.end(); p0++ ) |
|
|
819 // if (arg0.key (p0) == arg1) |
|
|
820 // @{ |
|
|
821 // tmp = arg0.contents (p0) (0); |
|
|
822 // break; |
|
|
823 // @} |
|
|
824 // though using seek is more concise. |
|
|
825 Octave_map::const_iterator p1 = arg0.seek (arg1); |
|
|
826 octave_value tmp = arg0.contents( p1 ) (0); |
|
|
827 Octave_map st; |
|
|
828 st.assign ("selected", tmp); |
|
|
829 retval = octave_value (st); |
|
|
830 @} |
|
|
831 @} |
|
|
832 return retval; |
|
|
833 @} |
|
|
834 @end group |
|
|
835 @end example |
|
|
836 |
|
|
837 An example of its use is |
|
|
838 |
|
|
839 @example |
|
|
840 @group |
|
|
841 x.a = 1; x.b = "test"; x.c = [1,2]; |
|
|
842 structdemo(x,"b") |
|
|
843 @result{} selected = test |
|
|
844 @end group |
|
|
845 @end example |
|
|
846 |
|
|
847 The commented code above demonstrates how to iterated over all of the |
|
|
848 fields of the structure, where as the following code demonstrates finding |
|
|
849 a particular field in a more concise manner. |
|
|
850 |
|
|
851 As can be seen the @code{contents} method of the @code{Octave_map} class |
|
|
852 returns a Cell which allows Structure Arrays to be |
|
|
853 represented. Therefore, to obtain the underlying octave_value we write |
|
|
854 |
|
|
855 @example |
|
|
856 octave_value tmp = arg0.contents (p1) (0); |
|
|
857 @end example |
|
|
858 |
|
|
859 where the trailing (0) is the () operator on the Cell array. |
|
|
860 |
|
|
861 @node Accessing Global Variables in Oct-Files |
|
|
862 @subsection Accessing Global Variables in Oct-Files |
|
|
863 |
|
|
864 Global variables allow variables in the global scope to be |
|
|
865 accessed. Global variables can easily be accessed with oct-files using |
|
|
866 the support functions @code{get_global_value} and |
|
|
867 @code{set_global_value}. @code{get_global_value} takes two arguments, |
|
|
868 the first is a string representing the variable name to obtain. The |
|
|
869 second argument is a boolean argument specifying what to do in the case |
|
|
870 that no global variable of the desired name is found. An example of the |
|
|
871 use of these two functions is |
|
|
872 |
|
|
873 @example |
|
|
874 @group |
|
|
875 #include <octave/oct.h> |
|
|
876 |
|
|
877 DEFUN_DLD (globaldemo, args, , "Global demo.") |
|
|
878 @{ |
|
|
879 int nargin = args.length(); |
|
|
880 octave_value retval; |
|
|
881 |
|
|
882 if (nargin != 1) |
|
|
883 print_usage (); |
|
|
884 else |
|
|
885 @{ |
|
|
886 std::string s = args(0).string_value (); |
|
|
887 if (! error_state) |
|
|
888 @{ |
|
|
889 octave_value tmp = get_global_value (s, true); |
|
|
890 if (tmp.is_defined ()) |
|
|
891 retval = tmp; |
|
|
892 else |
|
|
893 retval = "Global variable not found"; |
|
|
894 |
|
|
895 set_global_value ("a", 42.0); |
|
|
896 @} |
|
|
897 @} |
|
|
898 return retval; |
|
|
899 @} |
|
|
900 @end group |
|
|
901 @end example |
|
|
902 |
|
|
903 An example of its use is |
|
|
904 |
|
|
905 @example |
|
|
906 @group |
|
|
907 global a b |
|
|
908 b = 10; |
|
|
909 globaldemo ("b") |
|
|
910 @result{} 10 |
|
|
911 globaldemo ("c") |
|
|
912 @result{} "Global variable not found" |
|
|
913 num2str (a) |
|
|
914 @result{} 42 |
|
|
915 @end group |
|
|
916 @end example |
|
|
917 |
|
|
918 @node Calling Octave Functions from Oct-Files |
|
|
919 @subsection Calling Octave Functions from Oct-Files |
|
|
920 |
|
|
921 There is often a need to be able to call another octave function from |
|
|
922 within an oct-file, and there are many examples of such within octave |
|
|
923 itself. For example the @code{quad} function is an oct-file that |
|
|
924 calculates the definite integral by quadrature over a user supplied |
|
|
925 function. |
|
|
926 |
|
|
927 There are also many ways in which a function might be passed. It might |
|
|
928 be passed as one of |
|
|
929 |
|
|
930 @enumerate 1 |
|
|
931 @item Function Handle |
|
|
932 @item Anonymous Function Handle |
|
|
933 @item Inline Function |
|
|
934 @item String |
|
|
935 @end enumerate |
|
|
936 |
|
|
937 The example below demonstrates an example that accepts all four means of |
|
|
938 passing a function to an oct-file. |
|
|
939 |
|
|
940 @example |
|
|
941 @group |
|
|
942 #include <octave/oct.h> |
|
|
943 #include <octave/parse.h> |
|
|
944 |
|
|
945 DEFUN_DLD (funcdemo, args, nargout, "Function Demo") |
|
|
946 @{ |
|
|
947 int nargin = args.length(); |
|
|
948 octave_value_list retval; |
|
|
949 |
|
|
950 if (nargin < 2) |
|
|
951 print_usage (); |
|
|
952 else |
|
|
953 @{ |
|
|
954 octave_value_list newargs; |
|
|
955 for (octave_idx_type i = nargin - 1; i > 0; i--) |
|
|
956 newargs (i - 1) = args(i); |
|
|
957 if (args(0).is_function_handle () || |
|
|
958 args(0).is_inline_function ()) |
|
|
959 @{ |
|
|
960 octave_function *fcn = args(0).function_value (); |
|
|
961 if (! error_state) |
|
|
962 retval = feval (fcn, newargs, nargout); |
|
|
963 @} |
|
|
964 else if (args(0).is_string ()) |
|
|
965 @{ |
|
|
966 std::string fcn = args (0).string_value (); |
|
|
967 if (! error_state) |
|
|
968 retval = feval (fcn, newargs, nargout); |
|
|
969 @} |
|
|
970 else |
|
|
971 error ("funcdemo: expected string, inline or function handle"); |
|
|
972 @} |
|
|
973 return retval; |
|
|
974 @} |
|
|
975 @end group |
|
|
976 @end example |
|
|
977 |
|
|
978 The first argument to this demonstration is the user supplied function |
|
|
979 and the following arguments are all passed to the user function. |
|
|
980 |
|
|
981 @example |
|
|
982 @group |
|
|
983 funcdemo(@@sin,1) |
|
|
984 @result{} 0.84147 |
|
|
985 funcdemo(@@(x) sin(x), 1) |
|
|
986 @result{} 0.84147 |
|
|
987 funcdemo (inline("sin(x)"), 1) |
|
|
988 @result{} 0.84147 |
|
|
989 funcdemo("sin",1) |
|
|
990 @result{} 0.84147 |
|
|
991 funcdemo (@@atan2, 1, 1) |
|
|
992 @result{} 0.78540 |
|
|
993 @end group |
|
|
994 @end example |
|
|
995 |
|
|
996 When the user function is passed as a string, the treatment of the |
|
|
997 function is different. In some cases it is necessary to always have the |
|
|
998 user supplied function as an octave_function. In that case the string |
|
|
999 argument can be used to create a temporary function like |
|
|
1000 |
|
|
1001 @example |
|
|
1002 @group |
|
|
1003 std::octave fcn_name = unique_symbol_name ("__fcn__"); |
|
|
1004 std::string fname = "function y = "; |
|
|
1005 fname.append (fcn_name); |
|
|
1006 fname.append ("(x) y = "); |
|
|
1007 fcn = extract_function (args(0), "funcdemo", fcn_name, |
|
|
1008 fname, "; endfunction"); |
|
|
1009 @dots{} |
|
|
1010 if (fcn_name.length()) |
|
|
1011 clear_function (fcn_name); |
|
|
1012 @end group |
|
|
1013 @end example |
|
|
1014 |
|
|
1015 There are two important things to know in this case. The number of input |
|
|
1016 arguments to the user function is fixed, and in the above is a single |
|
|
1017 argument, and secondly to avoid leaving the temporary function in the |
|
|
1018 Octave symbol table it should be cleared after use. |
|
|
1019 |
|
|
1020 @node Calling External Code from Oct-Files |
|
|
1021 @subsection Calling External Code from Oct-Files |
|
|
1022 |
|
|
1023 Linking external C code to Octave is relatively simple, as the C |
|
|
1024 functions can easily be called directly from C++. One possible issue is |
|
|
1025 the declarations of the external C functions might need to be explicitly |
|
|
1026 defined as C functions to the compiler. If the declarations of the |
|
|
1027 external C functions are in the header @code{foo.h}, then the manner in |
|
|
1028 which to ensure that the C++ compiler treats these declarations as C |
|
|
1029 code is |
|
|
1030 |
|
|
1031 @example |
|
|
1032 @group |
|
|
1033 #ifdef __cplusplus |
|
|
1034 extern "C" |
|
|
1035 @{ |
|
|
1036 #endif |
|
|
1037 #include "foo.h" |
|
|
1038 #ifdef __cplusplus |
|
|
1039 @} /* end extern "C" */ |
|
|
1040 #endif |
|
|
1041 @end group |
|
|
1042 @end example |
|
|
1043 |
|
|
1044 Calling Fortran code however can pose some difficulties. This is due to |
|
|
1045 differences in the manner in compilers treat the linking of Fortran code |
|
|
1046 with C or C++ code. Octave supplies a number of macros that allow |
|
|
1047 consistent behavior across a number of compilers. |
|
|
1048 |
|
|
1049 The underlying Fortran code should use the @code{XSTOPX} function to |
|
|
1050 replace the Fortran @code{STOP} function. @code{XSTOPX} uses the Octave |
|
|
1051 exception handler to treat failing cases in the fortran code |
|
|
1052 explicitly. Note that Octave supplies its own replacement blas |
|
|
1053 @code{XERBLA} function, which uses @code{XSTOPX}. |
|
|
1054 |
|
|
1055 If the underlying code calls @code{XSTOP}, then the @code{F77_XFCN} |
|
|
1056 macro should be used to call the underlying fortran function. The Fortran |
|
|
1057 exception state can then be checked with the global variable |
|
|
1058 @code{f77_exception_encountered}. If @code{XSTOP} will not be called, |
|
|
1059 then the @code{F77_FCN} macro should be used instead to call the Fortran |
|
|
1060 code. |
|
|
1061 |
|
|
1062 There is no harm in using @code{F77_XFCN} in all cases, except that for |
|
|
1063 Fortran code that is short running and executes a large number of times, |
|
|
1064 there is potentially an overhead in doing so. However, if @code{F77_FCN} |
|
|
1065 is used with code that calls @code{XSTOP}, Octave can generate a |
|
|
1066 segmentation fault. |
|
|
1067 |
|
|
1068 An example of the inclusion of a Fortran function in an oct-file is |
|
|
1069 given in the following example, where the C++ wrapper is |
|
|
1070 |
|
|
1071 @example |
|
|
1072 @group |
|
|
1073 #include <octave/oct.h> |
|
|
1074 #include <octave/f77-fcn.h> |
|
|
1075 |
|
|
1076 extern "C" |
|
|
1077 @{ |
|
|
1078 F77_RET_T |
|
|
1079 F77_FUNC (fortsub, FORTSUB) |
|
|
1080 (const int&, double*, F77_CHAR_ARG_DECL |
|
|
1081 F77_CHAR_ARG_LEN_DECL); |
|
|
1082 @} |
|
|
1083 |
|
|
1084 DEFUN_DLD (fortdemo , args , , "Fortran Demo.") |
|
|
1085 @{ |
|
|
1086 octave_value_list retval; |
|
|
1087 int nargin = args.length(); |
|
|
1088 if (nargin != 1) |
|
|
1089 print_usage (); |
|
|
1090 else |
|
|
1091 @{ |
|
|
1092 NDArray a = args(0).array_value (); |
|
|
1093 if (! error_state) |
|
|
1094 @{ |
|
|
1095 double *av = a.fortran_vec (); |
|
|
1096 octave_idx_type na = a.nelem (); |
|
|
1097 OCTAVE_LOCAL_BUFFER (char, ctmp, 128); |
|
|
1098 |
|
|
1099 F77_XFCN(fortsub, FORTSUB, (na, av, ctmp F77_CHAR_ARG_LEN (128))); |
|
|
1100 |
|
|
1101 if ( f77_exception_encountered ) |
|
|
1102 error ("fortdemo: error in fortran"); |
|
|
1103 else |
|
|
1104 @{ |
|
|
1105 retval(1) = std::string (ctmp); |
|
|
1106 retval(0) = a; |
|
|
1107 @} |
|
|
1108 @} |
|
|
1109 @} |
|
|
1110 return retval; |
|
|
1111 @} |
|
|
1112 @end group |
|
|
1113 @end example |
|
|
1114 |
|
|
1115 and the fortran function is |
|
|
1116 |
|
|
1117 @example |
|
|
1118 @group |
|
|
1119 subroutine fortsub (n, a, s) |
|
|
1120 implicit none |
|
|
1121 character*(*) s |
|
|
1122 real*8 a(*) |
|
|
1123 integer*4 i, n, ioerr |
|
|
1124 do i = 1, n |
|
|
1125 if (a (i) .eq. 0d0) then |
|
|
1126 call xstopx('fortsub:divide by zero') |
|
|
1127 else |
|
|
1128 a (i) = 1d0 / a (i) |
|
|
1129 end if |
|
|
1130 end do |
|
|
1131 write(unit = s, fmt = '(a,i3,a,a)', iostat = ioerr) |
|
|
1132 1 'There are ', n, ' values in the input vector', char(0) |
|
|
1133 if (ioerr .ne. 0) then |
|
|
1134 call xstopx('fortsub: error writing string') |
|
|
1135 end if |
|
|
1136 return |
|
|
1137 end |
|
|
1138 @end group |
|
|
1139 @end example |
|
|
1140 |
|
|
1141 This example demonstrates most of the features needed to link to an |
|
|
1142 external Fortran function, including passing arrays and strings, as well |
|
|
1143 as exception handling. An example of the behavior of this function is |
|
|
1144 |
|
|
1145 @example |
|
|
1146 @group |
|
|
1147 [b, s] = fortdemo(1:3) |
|
|
1148 @result{} |
|
|
1149 b = 1.00000 0.50000 0.33333 |
|
|
1150 s = There are 3 values in the input vector |
|
|
1151 [b, s] = fortdemo(0:3) |
|
|
1152 error: fortsub:divide by zero |
|
|
1153 error: exception encountered in Fortran subroutine fortsub_ |
|
|
1154 error: fortdemo: error in fortran |
|
|
1155 @end group |
|
|
1156 @end example |
|
|
1157 |
|
|
1158 @node Allocating Local Memory in Oct-Files |
|
|
1159 @subsection Allocating Local Memory in Oct-Files |
|
|
1160 |
|
|
1161 Allocating memory within an oct-file might seem easy as the C++ |
|
|
1162 new/delete operators can be used. However, in that case care must be |
|
|
1163 taken to avoid memory leaks. The preferred manner in which to allocate |
|
|
1164 memory for use locally is to use the OCTAVE_LOCAL_BUFFER macro. An |
|
|
1165 example of its use is |
|
|
1166 |
|
|
1167 @example |
|
|
1168 OCTAVE_LOCAL_BUFFER (double, tmp, len) |
|
|
1169 @end example |
|
|
1170 |
|
|
1171 that returns a pointer @code{tmp} of type @code{double *} of length |
|
|
1172 @code{len}. |
|
|
1173 |
|
|
1174 @node Input Parameter Checking in Oct-Files |
|
|
1175 @subsection Input Parameter Checking in Oct-Files |
|
|
1176 |
|
|
1177 WRITE ME |
|
|
1178 |
|
|
1179 @node Exception and Error Handling in Oct-Files |
|
|
1180 @subsection Exception and Error Handling in Oct-Files |
|
|
1181 |
|
|
1182 Another important feature of Octave is its ability to react to the user |
|
|
1183 typing @kbd{Control-C} even during calculations. This ability is based on the |
|
|
1184 C++ exception handler, where memory allocated by the C++ new/delete |
|
|
1185 methods are automatically released when the exception is treated. When |
|
|
1186 writing an oct-file, to allow Octave to treat the user typing @kbd{Control-C}, |
|
|
1187 the @code{OCTAVE_QUIT} macro is supplied. For example |
|
|
1188 |
|
|
1189 @example |
|
|
1190 @group |
|
|
1191 for (octave_idx_type i = 0; i < a.nelem (); i++) |
|
|
1192 @{ |
|
|
1193 OCTAVE_QUIT; |
|
|
1194 b.elem(i) = 2. * a.elem(i); |
|
|
1195 @} |
|
|
1196 @end group |
|
|
1197 @end example |
|
|
1198 |
|
|
1199 The presence of the OCTAVE_QUIT macro in the inner loop allows Octave to |
|
|
1200 treat the user request with the @kbd{Control-C}. Without this macro, the user |
|
|
1201 must either wait for the function to return before the interrupt is |
|
|
1202 processed, or press @kbd{Control-C} three times to force Octave to exit. |
|
|
1203 |
|
|
1204 The OCTAVE_QUIT macro does impose a very small speed penalty, and so for |
|
|
1205 loops that are known to be small it might not make sense to include |
|
|
1206 OCTAVE_QUIT. |
|
|
1207 |
|
|
1208 When creating an oct-file that uses an external libraries, the function |
|
|
1209 might spend a significant portion of its time in the external |
|
|
1210 library. It is not generally possible to use the OCTAVE_QUIT macro in |
|
|
1211 this case. The alternative in this case is |
|
|
1212 |
|
|
1213 @example |
|
|
1214 @group |
|
|
1215 BEGIN_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE; |
|
|
1216 @dots{} some code that calls a "foreign" function @dots{} |
|
|
1217 END_INTERRUPT_IMMEDIATELY_IN_FOREIGN_CODE; |
|
|
1218 @end group |
|
|
1219 @end example |
|
|
1220 |
|
|
1221 The disadvantage of this is that if the foreign code allocates any |
|
|
1222 memory internally, then this memory might be lost during an interrupt, |
|
|
1223 without being deallocated. Therefore, ideally Octave itself should |
|
|
1224 allocate any memory that is needed by the foreign code, with either the |
|
|
1225 fortran_vec method or the OCTAVE_LOCAL_BUFFER macro. |
|
|
1226 |
|
|
1227 The Octave unwind_protect mechanism (@ref{The unwind_protect Statement}) |
|
|
1228 can also be used in oct-files. In conjunction with the exception |
|
|
1229 handling of Octave, it is important to enforce that certain code is run |
|
|
1230 to allow variables, etc to be restored even if an exception occurs. An |
|
|
1231 example of the use of this mechanism is |
|
|
1232 |
|
|
1233 @example |
|
|
1234 @group |
|
|
1235 #include <octave/oct.h> |
|
|
1236 #include <octave/unwind-prot.h> |
|
|
1237 |
|
|
1238 void |
|
|
1239 err_hand (const char *fmt, ...) |
|
|
1240 @{ |
|
|
1241 // Do nothing!! |
|
|
1242 @} |
|
|
1243 |
|
|
1244 DEFUN_DLD (unwinddemo, args, nargout, "Unwind Demo") |
|
|
1245 @{ |
|
|
1246 int nargin = args.length(); |
|
|
1247 octave_value retval; |
|
|
1248 if (nargin < 2) |
|
|
1249 print_usage (); |
|
|
1250 else |
|
|
1251 @{ |
|
|
1252 NDArray a = args(0).array_value (); |
|
|
1253 NDArray b = args(1).array_value (); |
|
|
1254 |
|
|
1255 if (! error_state) |
|
|
1256 @{ |
|
|
1257 unwind_protect::begin_frame("Funwinddemo"); |
|
|
1258 unwind_protect_ptr(current_liboctave_warning_handler); |
|
|
1259 set_liboctave_warning_handler(err_hand); |
|
|
1260 retval = octave_value ( quotient (a, b)); |
|
|
1261 unwind_protect::run_frame("Funwinddemo"); |
|
|
1262 @} |
|
|
1263 @} |
|
|
1264 return retval; |
|
|
1265 @} |
|
|
1266 @end group |
|
|
1267 @end example |
|
|
1268 |
|
|
1269 As can be seen in the example |
|
|
1270 |
|
|
1271 @example |
|
|
1272 @group |
|
|
1273 unwinddemo(1,0) |
|
|
1274 @result{} Inf |
|
|
1275 1 / 0 |
|
|
1276 @result{} warning: division by zero |
|
|
1277 Inf |
|
|
1278 @end group |
|
|
1279 @end example |
|
|
1280 |
|
|
1281 The division by zero (and in fact all warnings) is disabled in the |
|
|
1282 @code{unwinddemo} function. |
|
|
1283 |
|
|
1284 @node Documentation and Test of Oct-Files |
|
|
1285 @subsection Documentation and Test of Oct-Files |
|
|
1286 |
|
|
1287 WRITE ME, reference how to use texinfo in oct-file and embed test code. |
|
|
1288 |
|
|
1289 @node Application Programming Interface for Oct-Files |
|
|
1290 @subsection Application Programming Interface for Oct-Files |
|
|
1291 |
|
|
1292 WRITE ME, using Coda section 1.3 as a starting point. |
|
|
1293 |
|
|
1294 @node Mex-Files |
|
|
1295 @section Mex-Files |
|
|
1296 @cindex mex-files |
|
|
1297 @cindex mex |
|
|
1298 |
|
|
1299 Octave includes an interface to allow legacy mex-files to be compiled |
|
|
1300 and used with Octave. This interface can also be used to share code |
|
|
1301 between Octave and non Octave users. However, as mex-files expose the |
|
|
1302 intern API of a product alternative to Octave, and the internal |
|
|
1303 structure of Octave is different to this product, a mex-file can never |
|
|
1304 have the same performance in Octave as the equivalent oct-file. In |
|
|
1305 particular to support the manner in which mex-files access the variables |
|
|
1306 passed to mex functions, there are a significant number of additional |
|
|
1307 copies of memory when calling or returning from a mex function. For this |
|
|
1308 reason, new code should be written using the oct-file interface |
|
|
1309 discussed above if possible. |
|
|
1310 |
|
|
1311 @menu |
|
|
1312 * Getting Started with Mex-Files:: |
|
|
1313 * Using Structures with Mex-Files:: |
|
|
1314 * Sparse Matrices with Mex-Files:: |
|
|
1315 * Calling External Functions in Mex-Files:: |
|
|
1316 @end menu |
|
|
1317 |
|
|
1318 @node Getting Started with Mex-Files |
|
|
1319 @subsection Getting Started with Mex-Files |
|
|
1320 |
|
|
1321 The basic command to build a mex-file is either @code{mkoctfile --mex} or |
|
|
1322 @code{mex}. The first can either be used from within Octave or from the |
|
|
1323 commandline. However, to avoid issues with the installation of other |
|
|
1324 products, the use of the command @code{mex} is limited to within Octave. |
|
|
1325 |
|
|
1326 @DOCSTRING(mex) |
|
|
1327 |
|
|
1328 @DOCSTRING(mexext) |
|
|
1329 |
|
|
1330 One important difference between the use of mex with other products and |
|
|
1331 with Octave is that the header file "matrix.h" is implicitly included |
|
|
1332 through the inclusion of "mex.h". This is to avoid a conflict with the |
|
|
1333 Octave file "Matrix.h" with operating systems and compilers that don't |
|
|
1334 distinguish between filenames in upper and lower case |
|
|
1335 |
|
|
1336 Consider the short example |
|
|
1337 |
|
|
1338 @example |
|
|
1339 @group |
|
|
1340 #include "mex.h" |
|
|
1341 |
|
|
1342 void |
|
|
1343 mexFunction (int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) |
|
|
1344 @{ |
|
|
1345 mxArray *v = mxCreateDoubleMatrix (1, 1, mxREAL); |
|
|
1346 double *data = mxGetPr (v); |
|
|
1347 *data = 1.23456789; |
|
|
1348 plhs[0] = v; |
|
|
1349 @} |
|
|
1350 @end group |
|
|
1351 @end example |
|
|
1352 |
|
|
1353 This simple example demonstrates the basics of writing a mex-file. |
|
|
1354 |
|
|
1355 WRITE ME |
|
|
1356 |
|
|
1357 @node Using Structures with Mex-Files |
|
|
1358 @subsection Using Structures with Mex-Files |
|
|
1359 |
|
|
1360 WRITE ME |
|
|
1361 |
|
|
1362 @node Sparse Matrices with Mex-Files |
|
|
1363 @subsection Sparse Matrices with Mex-Files |
|
|
1364 |
|
|
1365 WRITE ME |
|
|
1366 |
|
|
1367 @node Calling External Functions in Mex-Files |
|
|
1368 @subsection Calling External Functions in Mex-Files |
|
|
1369 |
|
|
1370 WRITE ME |
|
|
1371 |
|
|
1372 @node Standalone Programs |
|
|
1373 @section Standalone Programs |
|
|
1374 |
|
|
1375 The libraries Octave itself uses, can be utilized in standalone |
|
|
1376 applications. These applications then have access, for example, to the |
|
|
1377 array and matrix classes as well as to all the Octave algorithms. The |
|
|
1378 following C++ program, uses class Matrix from liboctave.a or |
|
|
1379 liboctave.so. |
|
|
1380 |
|
|
1381 @example |
|
|
1382 @group |
|
|
1383 #include <iostream> |
|
|
1384 #include <octave/oct.h> |
|
|
1385 int |
|
|
1386 main (void) |
|
|
1387 @{ |
|
|
1388 std::cout << "Hello Octave world!\n"; |
|
|
1389 const int size = 2; |
|
|
1390 Matrix a_matrix = Matrix(size, size); |
|
|
1391 for (octave_idx_type row = 0; row < size; ++row) |
|
|
1392 @{ |
|
|
1393 for (octave_idx_type column = 0; column < size; ++column) |
|
|
1394 @{ |
|
|
1395 a_matrix(row, column) = (row + 1)*10 + (column + 1); |
|
|
1396 @} |
|
|
1397 @} |
|
|
1398 std::cout << a_matrix; |
|
|
1399 return 0; |
|
|
1400 @} |
|
|
1401 @end group |
|
|
1402 @end example |
|
|
1403 |
|
|
1404 mkoctfile can then be used to build a standalone application with a |
|
|
1405 command like |
|
|
1406 |
|
|
1407 @example |
|
|
1408 @group |
|
|
1409 $ mkoctfile --link-stand-alone hello.cc -o hello |
|
|
1410 $ ./hello |
|
|
1411 Hello Octave world! |
|
|
1412 11 12 |
|
|
1413 21 22 |
|
|
1414 $ |
|
|
1415 @end group |
|
|
1416 @end example |
|
|
1417 |
|
|
1418 Note that the application @code{hello} will be dynamically linked |
|
|
1419 against the octave libraries and any octave support libraries. |
