SCM

SCM Repository

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

Diff of /pkg/Matrix/src/sparseQR.c

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

revision 1432, Thu Aug 24 16:21:18 2006 UTC revision 1960, Fri Jul 6 16:54:43 2007 UTC
# Line 2  Line 2 
2    
3  SEXP sparseQR_validate(SEXP x)  SEXP sparseQR_validate(SEXP x)
4  {  {
5      cs *V = Matrix_as_cs(GET_SLOT(x, install("V"))),      CSP V = AS_CSP(GET_SLOT(x, install("V"))),
6          *R = Matrix_as_cs(GET_SLOT(x, install("R")));          R = AS_CSP(GET_SLOT(x, install("R")));
7      SEXP beta = GET_SLOT(x, install("beta")),      SEXP beta = GET_SLOT(x, install("beta")),
8          p = GET_SLOT(x, Matrix_pSym),          p = GET_SLOT(x, Matrix_pSym),
9          q = GET_SLOT(x, install("q"));          q = GET_SLOT(x, install("q"));
10      int lq = LENGTH(q);      int lq = LENGTH(q);
11        R_CheckStack();
12    
13      if (LENGTH(p) != V->m)      if (LENGTH(p) != V->m)
14          return mkString(_("length(p) must match nrow(V)"));          return mkString(_("length(p) must match nrow(V)"));
# Line 37  Line 38 
38                      double *y, int *ydims)                      double *y, int *ydims)
39  {  {
40      int j, k, m = V->m, n = V->n;      int j, k, m = V->m, n = V->n;
41      double *x = Calloc(m, double);      /* workspace */      double *x = Alloca(m, double);      /* workspace */
42        R_CheckStack();
43    
44      if (ydims[0] != m)      if (ydims[0] != m)
45          error(_("Dimensions of system are inconsistent"));          error(_("Dimensions of system are inconsistent"));
# Line 55  Line 57 
57              Memcpy(yj, x, m);              Memcpy(yj, x, m);
58          }          }
59      }      }
     Free(x);  
60  }  }
61    
62    
63  SEXP sparseQR_qty(SEXP qr, SEXP y, SEXP trans)  SEXP sparseQR_qty(SEXP qr, SEXP y, SEXP trans)
64  {  {
65      SEXP ans = PROTECT(dup_mMatrix_as_dgeMatrix(y));      SEXP ans = PROTECT(dup_mMatrix_as_dgeMatrix(y));
66      cs *V = Matrix_as_cs(GET_SLOT(qr, install("V")));      CSP V = AS_CSP(GET_SLOT(qr, install("V")));
67        R_CheckStack();
68    
69      sparseQR_Qmult(V, REAL(GET_SLOT(qr, install("beta"))),      sparseQR_Qmult(V, REAL(GET_SLOT(qr, install("beta"))),
70                     INTEGER(GET_SLOT(qr, Matrix_pSym)),                     INTEGER(GET_SLOT(qr, Matrix_pSym)),
71                     asLogical(trans),                     asLogical(trans),
72                     REAL(GET_SLOT(ans, Matrix_xSym)),                     REAL(GET_SLOT(ans, Matrix_xSym)),
73                     INTEGER(GET_SLOT(ans, Matrix_DimSym)));                     INTEGER(GET_SLOT(ans, Matrix_DimSym)));
   
     Free(V);  
74      UNPROTECT(1);      UNPROTECT(1);
75      return ans;      return ans;
76  }  }
# Line 79  Line 79 
79  {  {
80      SEXP ans = PROTECT(dup_mMatrix_as_dgeMatrix(y)),      SEXP ans = PROTECT(dup_mMatrix_as_dgeMatrix(y)),
81          qslot = GET_SLOT(qr, install("q"));          qslot = GET_SLOT(qr, install("q"));
82      cs *V = Matrix_as_cs(GET_SLOT(qr, install("V"))),      CSP V = AS_CSP(GET_SLOT(qr, install("V"))),
83          *R = Matrix_as_cs(GET_SLOT(qr, install("R")));          R = AS_CSP(GET_SLOT(qr, install("R")));
84      int *ydims = INTEGER(GET_SLOT(ans, Matrix_DimSym)),      int *ydims = INTEGER(GET_SLOT(ans, Matrix_DimSym)),
85          *q = INTEGER(qslot),          *q = INTEGER(qslot),
86          j, lq = LENGTH(qslot), m = R->m, n = R->n;          j, lq = LENGTH(qslot), m = R->m, n = R->n;
87      double *ax = REAL(GET_SLOT(ans, Matrix_xSym)),      double *ax = REAL(GET_SLOT(ans, Matrix_xSym)),
88          *x = Calloc(m, double);          *x = Alloca(m, double);
89        R_CheckStack();
90        R_CheckStack();
91    
92      /* apply row permutation and multiply by Q' */      /* apply row permutation and multiply by Q' */
93      sparseQR_Qmult(V, REAL(GET_SLOT(qr, install("beta"))),      sparseQR_Qmult(V, REAL(GET_SLOT(qr, install("beta"))),
# Line 99  Line 101 
101              Memcpy(aj, x, n);              Memcpy(aj, x, n);
102          }          }
103      }      }
     Free(V); Free(R); Free(x);  
104      UNPROTECT(1);      UNPROTECT(1);
105      return ans;      return ans;
106  }  }
# Line 107  Line 108 
108  SEXP sparseQR_resid_fitted(SEXP qr, SEXP y, SEXP resid)  SEXP sparseQR_resid_fitted(SEXP qr, SEXP y, SEXP resid)
109  {  {
110      SEXP ans = PROTECT(dup_mMatrix_as_dgeMatrix(y));      SEXP ans = PROTECT(dup_mMatrix_as_dgeMatrix(y));
111      cs *V = Matrix_as_cs(GET_SLOT(qr, install("V")));      CSP V = AS_CSP(GET_SLOT(qr, install("V")));
112      int *ydims = INTEGER(GET_SLOT(ans, Matrix_DimSym)),      int *ydims = INTEGER(GET_SLOT(ans, Matrix_DimSym)),
113          *p = INTEGER(GET_SLOT(qr, Matrix_pSym)),          *p = INTEGER(GET_SLOT(qr, Matrix_pSym)),
114          i, j, m = V->m, n = V->n, res = asLogical(resid);          i, j, m = V->m, n = V->n, res = asLogical(resid);
115      double *ax = REAL(GET_SLOT(ans, Matrix_xSym)),      double *ax = REAL(GET_SLOT(ans, Matrix_xSym)),
116          *beta = REAL(GET_SLOT(qr, install("beta")));          *beta = REAL(GET_SLOT(qr, install("beta")));
117        R_CheckStack();
118    
119      /* apply row permutation and multiply by Q' */      /* apply row permutation and multiply by Q' */
120      sparseQR_Qmult(V, beta, p, 1, ax, ydims);      sparseQR_Qmult(V, beta, p, 1, ax, ydims);
# Line 124  Line 126 
126      }      }
127      /* multiply by Q and apply inverse row permutation */      /* multiply by Q and apply inverse row permutation */
128      sparseQR_Qmult(V, beta, p, 0, ax, ydims);      sparseQR_Qmult(V, beta, p, 0, ax, ydims);
129    
130      UNPROTECT(1);      UNPROTECT(1);
131      return ans;      return ans;
132  }  }

Legend:
Removed from v.1432  
changed lines
  Added in v.1960

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