CODEWITHSUNDEEP
c++matlabfloating-pointprecisionmex
kiriloff
kiriloff

Reputation: 26335

Mex function called from Matlab and numerical difference

I have a Matlab scientific code, and I want to compare precisely the output from Matlab code containing a call to inversion function pinv() (native Matlab) to output of the same Matlab code calling now a C++ implementation of pinv().

I instrument the Matlab code to call the C++ pinv instead of the Matlab pinv(): I compile the C++ pinv to a mex file, using the relevant API.

Can the call to a function compiled to mex file introduce some numerical difference on top of the difference in implementation between Matlab and C++ pinv() ?

Upvotes: 0

Views: 250

Answers (1)

user2271770
user2271770

Reputation:

All the matrix operations executed by MATLAB are implemented by the LAPACK library so, yes, there might be differences induced by implementation. If what you're asking is "does the LAPACK library use different hardware floating point features than the C++ compiled code", that is highly unlikely—unless the C++ compiler is instructed to generate code targeted to a specific HW platform that subsets the available one.

However, even if the codes are similar at source level, LAPACK is compiled from FORTRAN source code, while the custom routine results from C++; the optimizations done by different compilers are different. So there's a high chance that - even with strikingly similar source code - the resulting executable code will be different. Which may lead to different results.

Note: The manipulation of data (i.e. the way is passed to a function) doesn't alter the content of the data; it may be a matter of execution time, not of numeric deviation. MEX is basically a convention on how the native MATLAB data structures are presented to whatever routine processes them; the efficiency of the routine accessing the data is again a matter of implementation.

Further comment: One must say that FORTRAN in-memory layout of 2D data is column-first, while the C-and-family in-memory layout of 2D data is row first (I'm talking about indexing). This can also induce differences in execution time (as in exec code also), depending on how the implementation algorithm accesses the data, because modern HW architectures' performance is highly dependent on, e.g. how many cache misses are occurring.

Upvotes: 1

Related Questions