[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: Cholesky updating
From: |
Jaroslav Hajek |
Subject: |
Re: Cholesky updating |
Date: |
Wed, 5 Mar 2008 15:32:54 +0100 |
Hello,
I've tried to address all the issues pointed out by John. A changeset
is attached.
Consider also applying the subsequent changeset - it changes the QR
row/col insert/delete methods to use 0-based indexing (for consistency
with the rest of liboctave, though this means slightly more fiddling
with indices).
On Wed, Mar 5, 2008 at 4:26 AM, John W. Eaton <address@hidden> wrote:
> On 2-Mar-2008, Jaroslav Hajek wrote:
>
> | please ignore the cholupdate changeset from the last message - it's
> | missing documentation and tests for cholupdate. A replacement is
> | attached.
> |
> | On Sat, Mar 1, 2008 at 9:42 PM, Jaroslav Hajek <address@hidden> wrote:
> | > hello,
> | >
> | > the attached changeset implements rank-1 updating for Cholesky
> | > factorization classes (CHOL, ComplexCHOL) in liboctave, as well as
> | > intrinsic function cholupdate, providing the functionality of the
> | > corresponding Matlab routine. The parent QR-updating changeset I've
> | > sent earlier is also attached for convenience.
>
> I applied this patch with some additional changes.
>
> instad of putting cholupdate in a separate file, I added it to
> src/DLD-FUNCTIONS/chol.cc.
>
> | + F77_RET_T
> | + F77_FUNC (dch1dn, DCH1DN) (const octave_idx_type&, double*, double*,
> double*,
> | + int &);
>
> All Fortran INTEGER values should be declared octave_idx_type so that
> --enable-64 can have a chance of working.
>
> | +void CHOL::update (Matrix u)
> | +{
> | + int n = chol_mat.rows ();
> | + if (u.length () != n)
> | + {
> | + (*current_liboctave_error_handler) ("CHOL update dimension
> mismatch");
> | + return;
> | + }
> | +
> | + OCTAVE_LOCAL_BUFFER(double, w, n);
> | +
> | + F77_XFCN(dch1up, DCH1UP, (n, chol_mat.fortran_vec (),
> | + u.fortran_vec (), w));
> | +}
>
> Instead of passing U by value here, I would prefer to pass it as
> const reference and then use an explicit temporary inside the
> function. Also, N should be declared octave_idx_type.
>
> | + octave_value argR,argu,argop;
> | +
> | + if ((nargin == 3 || nargin == 2) && nargout <= 2
> | + && (argR = args (0),argR.is_matrix_type ())
> | + && (argu = args (1),argu.is_matrix_type ())
> | + && (nargin < 3 || (argop = args(2), argop.is_string ())))
> | + {
> | + int n = argR.rows ();
> | + std::string op = (nargin < 3) ? "+" : argop.string_value();
> | + bool down = false;
> | +
> | + if (nargin < 3 || (op == "+") || (down = op == "-"))
> | + if (argR.columns () == n
> | + && argu.rows () == n && argu.columns () == 1)
> | + {
>
> Will you please fix this code to avoid comma expressions and
> assignments inside the if statement conditions?
>
> Your patch and my changes are now in my public archive. Will you
> please update and take a look at the argument parsing?
>
> Thanks,
>
> jwe
>
--
RNDr. Jaroslav Hajek
computing expert
Aeronautical Research and Test Institute (VZLU)
Prague, Czech Republic
url: www.highegg.matfyz.cz
qrchfixes.diff
Description: Text document
qr0based.diff
Description: Text document