On Tue, May 25, 2010 at 6:00 PM, Carlo de Falco <[email protected]> wrote: > > On 25 May 2010, at 16:17, Jaroslav Hajek wrote: > >> On Tue, May 25, 2010 at 3:13 PM, Carlo de Falco <[email protected]> >> wrote: >>> >>> On 25 May 2010, at 14:59, Jaroslav Hajek wrote: >>> >>>> unless anyone objects, I'll become the maintainer of this package. >>> >>> oh, well.. as now linear algebra has a maintainer I finally have someone >>> to >>> direct my complaints to ;) >>> >>>> comments, further suggestions? >>> >>> well yes, some time ago an initial implementation of pgmres was added to >>> the >>> package but the author[*] >>> seems to not have taken too much care for making it compatibile with >>> matlab's gmres function, maybe you >>> could have a look into this before the release? >>> >>> c. >>> >>> [*] OK, you got me, that lazy author is me... pgmres needs a lot of >>> improvement but it has been on my todo list >>> for a long time now, so although I think it would make sense to complete >>> it >>> before the release I will not at all >>> ask to delay the release because of this. >>> >> >> >> Is there particular reason why pgmres is a C++ file? > > the actual reason why pgmres was written in the first place was to > demonstrate the overloading of the '*' > operator in c++ and the use of a function instead of a matrix in the > matrix-free implementation of an iterative solver. > In principle I have no objections to switching to an m file if performance > is not affected, though. > >> The only inner >> loop I can see is the modified Gram-Schmidt orthogonalization: >> >> H.resize (reit+1, reit, 0.0); >> for (octave_idx_type j = 0; j < reit; j++) >> { >> H(j, reit-1) = (V.column (j).transpose ()) * tmp; >> tmp = tmp - (H (j, reit-1) * V.column (j)); >> } >> >> which would be useful to have as a separate function, such as: >> >> >> function [y, r] = mgorth (v, x); >> # orthogonalizes a column vector x relative to orthogonal basis v >> using the modified gram-schmidt orthogonalization. >> # mathematically, this is equivalent to >> # r = v'*x; y = x - v*r; >> # but computed in a more stable manner that avoids accumulating >> roundoff errors when done sequentially. >> >> >> pgmres could then be written as a m-file. What do you think? > > why is only the inner loop important, if many iterations are performed > doesn't the outer loop also impact speed? >
If a loop contains "heavy" operations, such as large matrix multiplication, matrix division, factorization, convolution, fft, etc., the interpreter overhead will typically become negligible compared to the time taken by these operations, so it doesn't matter whether it's written in m-code or C++. Also, when a callback function is called in the loop, it will most likely be an m-function, so the loop will involve interpreter overhead anyway. -- RNDr. Jaroslav Hajek, PhD computing expert & GNU Octave developer Aeronautical Research and Test Institute (VZLU) Prague, Czech Republic url: www.highegg.matfyz.cz ------------------------------------------------------------------------------ _______________________________________________ Octave-dev mailing list [email protected] https://lists.sourceforge.net/lists/listinfo/octave-dev
