SCM

SCM Repository

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

Diff of /pkg/src/geMatrix.c

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

revision 331, Fri Nov 12 21:03:24 2004 UTC revision 332, Fri Nov 12 21:04:36 2004 UTC
# Line 3  Line 3 
3  SEXP geMatrix_validate(SEXP obj)  SEXP geMatrix_validate(SEXP obj)
4  {  {
5      SEXP x = GET_SLOT(obj, Matrix_xSym),      SEXP x = GET_SLOT(obj, Matrix_xSym),
6          Dim = GET_SLOT(obj, Matrix_DimSym);          Dim = GET_SLOT(obj, Matrix_DimSym),
7            fact = GET_SLOT(obj, Matrix_factorization),
8            rc = GET_SLOT(obj, Matrix_rcondSym);
9      int m, n;      int m, n;
10    
11      if (length(Dim) != 2)      if (length(Dim) != 2)
12          return ScalarString(mkChar("Dim slot must have length 2"));          return ScalarString(mkChar("Dim slot must have length 2"));
13      m = INTEGER(Dim)[0]; n = INTEGER(Dim)[1];      m = INTEGER(Dim)[0]; n = INTEGER(Dim)[1];
14      if (m < 0 || n < 0)      if (m < 0 || n < 0)
15          return          return ScalarString(mkChar("Negative value(s) in Dim"));
             ScalarString(mkChar("Negative value(s) in Dim"));  
16      if (length(x) != m * n)      if (length(x) != m * n)
17          return          return ScalarString(mkChar("length of x slot != prod(Dim)"));
18              ScalarString(mkChar("x slot must have length that is prod(Dim)"));      if (length(fact) > 0 && getAttrib(fact, R_NamesSymbol) == R_NilValue)
19            return ScalarString(mkChar("factorization slot must be named list"));
20        if (length(rc) > 0 && getAttrib(rc, R_NamesSymbol) == R_NilValue)
21            return ScalarString(mkChar("rcond slot must be named numeric vector"));
22      return ScalarLogical(1);      return ScalarLogical(1);
23  }  }
24    
# Line 43  Line 47 
47  double set_rcond(SEXP obj, char *typstr)  double set_rcond(SEXP obj, char *typstr)
48  {  {
49      char typnm[] = {'\0', '\0'};      char typnm[] = {'\0', '\0'};
50      SEXP rcv = GET_SLOT(obj, install("rcond"));      SEXP rcv = GET_SLOT(obj, Matrix_rcondSym);
51      double rcond = get_double_by_name(rcv, typnm);      double rcond = get_double_by_name(rcv, typnm);
52    
53      typnm[0] = rcond_type(typstr);      typnm[0] = rcond_type(typstr);
# Line 59  Line 63 
63                           dims, &anorm, &rcond,                           dims, &anorm, &rcond,
64                           (double *) R_alloc(4*dims[0], sizeof(double)),                           (double *) R_alloc(4*dims[0], sizeof(double)),
65                           (int *) R_alloc(dims[0], sizeof(int)), &info);                           (int *) R_alloc(dims[0], sizeof(int)), &info);
66          SET_SLOT(obj, install("rcond"),          SET_SLOT(obj, Matrix_rcondSym,
67                   set_double_by_name(rcv, rcond, typnm));                   set_double_by_name(rcv, rcond, typnm));
68      }      }
69      return rcond;      return rcond;
# Line 75  Line 79 
79      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("poMatrix")));      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("poMatrix")));
80      int *Dims = INTEGER(GET_SLOT(x, Matrix_DimSym)),      int *Dims = INTEGER(GET_SLOT(x, Matrix_DimSym)),
81          *vDims = INTEGER(GET_SLOT(val, Matrix_DimSym));          *vDims = INTEGER(GET_SLOT(val, Matrix_DimSym));
82      int n = Dims[1];      int k = Dims[0], n = Dims[1];
83      double one = 1.0, zero = 0.0;      double one = 1.0, zero = 0.0;
84    
85      if ((*Dims) > 0 && n > 0) {      SET_SLOT(val, Matrix_factorization, allocVector(VECSXP, 0));
86        SET_SLOT(val, Matrix_rcondSym, allocVector(REALSXP, 0));
87        SET_SLOT(val, Matrix_uploSym, ScalarString(mkChar("U")));
88          vDims[0] = vDims[1] = n;          vDims[0] = vDims[1] = n;
89          SET_SLOT(val, Matrix_xSym, allocVector(REALSXP, n * n));          SET_SLOT(val, Matrix_xSym, allocVector(REALSXP, n * n));
90          F77_CALL(dsyrk)(CHAR(asChar(GET_SLOT(val, Matrix_uploSym))),      if (k > 0) {
91                         "T", Dims + 1, Dims, &one,          F77_CALL(dsyrk)("U", "T", &n, &k,
92                         REAL(GET_SLOT(x, Matrix_xSym)),                          &one, REAL(GET_SLOT(x, Matrix_xSym)), &k,
93                         Dims, &zero,                          &zero, REAL(GET_SLOT(val, Matrix_xSym)), &n);
94                         REAL(GET_SLOT(val, Matrix_xSym)), &n);      } else {
95            int i, nsqr = n * n;
96            double *xvals = REAL(GET_SLOT(val, Matrix_xSym));
97            for (i = 0; i < nsqr; i++) xvals[i] = 0.;
98      }      }
99    
100      UNPROTECT(1);      UNPROTECT(1);
101      return val;      return val;
102  }  }
# Line 96  Line 106 
106      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("geMatrix")));      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("geMatrix")));
107      int *xDims = INTEGER(GET_SLOT(x, Matrix_DimSym)),      int *xDims = INTEGER(GET_SLOT(x, Matrix_DimSym)),
108          *yDims = INTEGER(GET_SLOT(y, Matrix_DimSym)),          *yDims = INTEGER(GET_SLOT(y, Matrix_DimSym)),
109          *vDims = INTEGER(GET_SLOT(val, Matrix_DimSym));          *vDims;
110      int m = xDims[1], n = yDims[1];      int m = xDims[1], n = yDims[1];
111      double one = 1.0, zero = 0.0;      double one = 1.0, zero = 0.0;
112    
113        SET_SLOT(val, Matrix_rcondSym, allocVector(REALSXP, 0));
114        SET_SLOT(val, Matrix_factorization, allocVector(VECSXP, 0));
115        SET_SLOT(val, Matrix_DimSym, allocVector(INTSXP, 2));
116        vDims = INTEGER(GET_SLOT(val, Matrix_DimSym));
117      if ((*xDims) > 0 && (*yDims) > 0 && n > 0 && m > 0) {      if ((*xDims) > 0 && (*yDims) > 0 && n > 0 && m > 0) {
118          if (*xDims != *yDims)          if (*xDims != *yDims)
119              error("Dimensions of x and y are not compatible for crossprod");              error("Dimensions of x and y are not compatible for crossprod");
# Line 120  Line 134 
134      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("geMatrix")));      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("geMatrix")));
135      int *xDims = INTEGER(GET_SLOT(x, Matrix_DimSym)),      int *xDims = INTEGER(GET_SLOT(x, Matrix_DimSym)),
136          *yDims = INTEGER(getAttrib(y, R_DimSymbol)),          *yDims = INTEGER(getAttrib(y, R_DimSymbol)),
137          *vDims = INTEGER(GET_SLOT(val, Matrix_DimSym));          *vDims;
138      int m = xDims[1], n = yDims[1];      int m = xDims[1], n = yDims[1];
139      double one = 1.0, zero = 0.0;      double one = 1.0, zero = 0.0;
140    
141      if (!(isMatrix(y) && isReal(y)))      if (!(isMatrix(y) && isReal(y)))
142          error("Argument y must be a numeric (real) matrix");          error("Argument y must be a numeric (real) matrix");
143        SET_SLOT(val, Matrix_rcondSym, allocVector(REALSXP, 0));
144        SET_SLOT(val, Matrix_factorization, allocVector(VECSXP, 0));
145        SET_SLOT(val, Matrix_DimSym, allocVector(INTSXP, 2));
146        vDims = INTEGER(GET_SLOT(val, Matrix_DimSym));
147      if ((*xDims) > 0 && (*yDims) > 0 && n > 0 && m > 0) {      if ((*xDims) > 0 && (*yDims) > 0 && n > 0 && m > 0) {
148          if (*xDims != *yDims)          if (*xDims != *yDims)
149              error("Dimensions of x and y are not compatible for crossprod");              error("Dimensions of x and y are not compatible for crossprod");
# Line 219  Line 237 
237    
238    
239      if (dims[0] != dims[1]) error("Solve requires a square matrix");      if (dims[0] != dims[1]) error("Solve requires a square matrix");
240        SET_SLOT(val, Matrix_rcondSym, allocVector(REALSXP, 0));
241        SET_SLOT(val, Matrix_factorization, allocVector(VECSXP, 0));
242      SET_SLOT(val, Matrix_xSym, duplicate(GET_SLOT(lu, Matrix_xSym)));      SET_SLOT(val, Matrix_xSym, duplicate(GET_SLOT(lu, Matrix_xSym)));
243      x = REAL(GET_SLOT(val, Matrix_xSym));      x = REAL(GET_SLOT(val, Matrix_xSym));
244      SET_SLOT(val, Matrix_DimSym, duplicate(GET_SLOT(lu, Matrix_DimSym)));      SET_SLOT(val, Matrix_DimSym, duplicate(GET_SLOT(lu, Matrix_DimSym)));
# Line 245  Line 265 
265          error("Matrices are not conformable for multiplication");          error("Matrices are not conformable for multiplication");
266      if (m < 1 || n < 1 || k < 1)      if (m < 1 || n < 1 || k < 1)
267          error("Matrices with zero extents cannot be multiplied");          error("Matrices with zero extents cannot be multiplied");
268        SET_SLOT(val, Matrix_rcondSym, allocVector(REALSXP, 0));
269        SET_SLOT(val, Matrix_factorization, allocVector(VECSXP, 0));
270      SET_SLOT(val, Matrix_xSym, allocVector(REALSXP, m * n));      SET_SLOT(val, Matrix_xSym, allocVector(REALSXP, m * n));
271      SET_SLOT(val, Matrix_DimSym, allocVector(INTSXP, 2));      SET_SLOT(val, Matrix_DimSym, allocVector(INTSXP, 2));
272      cdims = INTEGER(GET_SLOT(val, Matrix_DimSym));      cdims = INTEGER(GET_SLOT(val, Matrix_DimSym));

Legend:
Removed from v.331  
changed lines
  Added in v.332

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