SCM

SCM Repository

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

Diff of /pkg/src/bCrosstab.c

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

revision 739, Tue May 17 17:42:04 2005 UTC revision 740, Tue May 17 17:42:35 2005 UTC
# Line 1  Line 1 
1  #include "bCrosstab.h"  #include "bCrosstab.h"
2  /* TODO:  /* TODO:
3   * - Only do a fill-reducing permutation on the first non-nested factor   * - Only do a fill-reducing permutation on the first non-nested factor?
  * - Alternatively: change the algorithm for the fill-reducing  
  *   permutation to a greedy or a picky algorithm.  
4   */   */
5    
6    
# Line 46  Line 44 
44      }      }
45      if (nextra) {      if (nextra) {
46          int cnr, ntot = cnz + nextra, pos;          int cnr, ntot = cnz + nextra, pos;
47          int          int *dims = INTEGER(getAttrib(GET_SLOT(C, Matrix_xSym), R_DimSymbol)),
             *dims = INTEGER(getAttrib(GET_SLOT(C, Matrix_xSym), R_DimSymbol)),  
48              *Ti = Memcpy((int *) Calloc(ntot, int), ci, cnz),              *Ti = Memcpy((int *) Calloc(ntot, int), ci, cnz),
49              *Tj = expand_cmprPt(anc, cp, Calloc(ntot, int)),              *Tj = expand_cmprPt(anc, cp, Calloc(ntot, int)),
50              *Ci = Calloc(ntot, int);              *Ci = Calloc(ntot, int);
# Line 83  Line 80 
80   * @param db pointer to the diagonal block   * @param db pointer to the diagonal block
81   * @param odb pointer to the off-diagonal block   * @param odb pointer to the off-diagonal block
82   */   */
83  static void  static R_INLINE void
84  diag_update(SEXP db, SEXP odb)  diag_update(SEXP db, SEXP odb, int n, int k)
85  {  {
     SEXP dbpP = GET_SLOT(db, Matrix_pSym),  
         odbpP = GET_SLOT(odb, Matrix_pSym);  
   
86      SET_SLOT(db, Matrix_iSym,      SET_SLOT(db, Matrix_iSym,
87               Matrix_lgCsyrk(1, 0, length(dbpP) - 1, length(odbpP) - 1,               Matrix_lgCsyrk(1, 0, n, k,
88                              INTEGER(GET_SLOT(odb, Matrix_iSym)),                              INTEGER(GET_SLOT(odb, Matrix_iSym)),
89                              INTEGER(odbpP), 1,                              INTEGER(GET_SLOT(odb, Matrix_pSym)),
90                                1,
91                              GET_SLOT(db, Matrix_iSym),                              GET_SLOT(db, Matrix_iSym),
92                              INTEGER(dbpP)));                              INTEGER(GET_SLOT(db, Matrix_pSym))));
93  }  }
94    
95  /**  /**
96   * Update an off-diagonal block of L from the blocked crosstabulation   * Update an off-diagonal block of L from the blocked crosstabulation
97   *   *
98   * @param A upper block   * @param A lower block
99   * @param B lower block   * @param B upper block
100   * @param C product block   * @param C product block
101   * @param nrA number of rows in A   * @param nrA number of rows in A
102   */   */
103  static void  static R_INLINE void
104  offdiag_update(SEXP A, SEXP B, SEXP C, int nrA)  offdiag_update(SEXP A, SEXP B, SEXP C, int m, int n, int k)
105  {  {
     SEXP ApP = GET_SLOT(A, Matrix_pSym),  
         CpP = GET_SLOT(C, Matrix_pSym);  
   
106      SET_SLOT(C, Matrix_iSym,      SET_SLOT(C, Matrix_iSym,
107               Matrix_lgClgCmm(0, 1, nrA, length(CpP) - 1,               Matrix_lgClgCmm(0, 1, m, n, k,
                              length(ApP) - 1,  
108                               INTEGER(GET_SLOT(A, Matrix_iSym)),                               INTEGER(GET_SLOT(A, Matrix_iSym)),
109                               INTEGER(ApP),                               INTEGER(GET_SLOT(A, Matrix_pSym)),
110                               INTEGER(GET_SLOT(B, Matrix_iSym)),                               INTEGER(GET_SLOT(B, Matrix_iSym)),
111                               INTEGER(GET_SLOT(B, Matrix_pSym)),                               INTEGER(GET_SLOT(B, Matrix_pSym)),
112                               1,                               1,
113                               GET_SLOT(C, Matrix_iSym),                               GET_SLOT(C, Matrix_iSym),
114                               INTEGER(CpP)));                               INTEGER(GET_SLOT(C, Matrix_pSym))));
115  }  }
116    
117          /* check dimensions of x slots and modify if necessary */          /* check dimensions of x slots and modify if necessary */
# Line 342  Line 333 
333                                           VECTOR_ELT(L, Lind(k,j)));                                           VECTOR_ELT(L, Lind(k,j)));
334          }          }
335          for (k = j+1; k < nf; k++) { /* Update remaining columns */          for (k = j+1; k < nf; k++) { /* Update remaining columns */
336              diag_update(VECTOR_ELT(ZZpO, k), VECTOR_ELT(L, Lind(k, j)));              diag_update(VECTOR_ELT(ZZpO, k), VECTOR_ELT(L, Lind(k, j)),
337                            nlev[k], nlev[j]);
338              for (i = k + 1; i < nf; i++)              for (i = k + 1; i < nf; i++)
339                  offdiag_update(VECTOR_ELT(L, Lind(k, j)),                  offdiag_update(VECTOR_ELT(L, Lind(i, j)),
340                                 VECTOR_ELT(L, Lind(i, j)),                                 VECTOR_ELT(L, Lind(k, j)),
341                                 VECTOR_ELT(L, Lind(i, k)),                                 VECTOR_ELT(L, Lind(i, k)),
342                                 nlev[k]);                                 nlev[i], nlev[k], nlev[j]);
343          }          }
344          check_x_slot_dim(VECTOR_ELT(ZZpO, j));          check_x_slot_dim(VECTOR_ELT(ZZpO, j));
345      }      }

Legend:
Removed from v.739  
changed lines
  Added in v.740

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