SCM

SCM Repository

[matrix] Diff of /pkg/src/dgeMatrix.c
ViewVC logotype

Diff of /pkg/src/dgeMatrix.c

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 479, Wed Feb 2 14:52:26 2005 UTC revision 495, Thu Feb 3 18:48:31 2005 UTC
# Line 465  Line 465 
465      UNPROTECT(1);      UNPROTECT(1);
466      return val;      return val;
467  }  }
468    
469    SEXP dgeMatrix_Schur(SEXP x, SEXP vectors)
470    {
471        int *dims = INTEGER(GET_SLOT(x, Matrix_DimSym));
472        int vecs = asLogical(vectors), info, izero = 0, lwork = -1, n = dims[0];
473        double *work, tmp;
474        char *nms[] = {"WR", "WI", "T", "Z", ""};
475        SEXP val = PROTECT(Matrix_make_named(VECSXP, nms));
476    
477        if (n != dims[1] || n < 1)
478            error("dgeMatrix_Schur: argument x must be a non-null square matrix");
479        SET_VECTOR_ELT(val, 0, allocVector(REALSXP, n));
480        SET_VECTOR_ELT(val, 1, allocVector(REALSXP, n));
481        SET_VECTOR_ELT(val, 2, allocMatrix(REALSXP, n, n));
482        Memcpy(REAL(VECTOR_ELT(val, 2)), REAL(GET_SLOT(x, Matrix_xSym)), n * n);
483        SET_VECTOR_ELT(val, 3, allocMatrix(REALSXP, vecs ? n : 0, vecs ? n : 0));
484        F77_CALL(dgees)(vecs ? "V" : "N", "N", NULL, dims, (double *) NULL, dims, &izero,
485                        (double *) NULL, (double *) NULL, (double *) NULL, dims,
486                        &tmp, &lwork, (int *) NULL, &info);
487        if (info) error("dgeMatrix_Schur: first call to dgees failed");
488        lwork = (int) tmp;
489        work = Calloc(lwork, double);
490        F77_CALL(dgees)(vecs ? "V" : "N", "N", NULL, dims, REAL(VECTOR_ELT(val, 2)), dims,
491                        &izero, REAL(VECTOR_ELT(val, 0)), REAL(VECTOR_ELT(val, 1)),
492                        REAL(VECTOR_ELT(val, 3)), dims, work, &lwork,
493                        (int *) NULL, &info);
494        if (info) error("dgeMatrix_Schur: dgees returned code %d", info);
495        Free(work);
496        UNPROTECT(1);
497        return val;
498    }

Legend:
Removed from v.479  
changed lines
  Added in v.495

root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge