test4test6 test7test8test5 R-Forge: Matrix Methods and Classes: SCM Repository
SCM

SCM Repository

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

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

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

revision 2111, Thu Jan 31 17:06:57 2008 UTC revision 2112, Mon Feb 18 08:24:46 2008 UTC
# Line 281  Line 281 
281   * Class of the value can be dtrMatrix or dgeMatrix   * Class of the value can be dtrMatrix or dgeMatrix
282   */   */
283  {  {
284        int k1 = asInteger(k1P), k2 = asInteger(k2P);
285    
286        if (k1 > k2) {
287            error(_("Lower band %d > upper band %d"), k1, k2);
288        }
289        else {
290      SEXP aa, ans = PROTECT(dup_mMatrix_as_dgeMatrix(x));      SEXP aa, ans = PROTECT(dup_mMatrix_as_dgeMatrix(x));
291      int *adims = INTEGER(GET_SLOT(ans, Matrix_DimSym)),      int *adims = INTEGER(GET_SLOT(ans, Matrix_DimSym)),
292          i, j, k1 = asInteger(k1P), k2 = asInteger(k2P);              i, j, m = adims[0], n = adims[1],
293      int m = adims[0], n = adims[1], sqr = (adims[0] == adims[1]),              sqr = (adims[0] == adims[1]),
294          tru = (k1 >= 0), trl = (k2 <= 0);          tru = (k1 >= 0), trl = (k2 <= 0);
295      double *ax = REAL(GET_SLOT(ans, Matrix_xSym));      double *ax = REAL(GET_SLOT(ans, Matrix_xSym));
296    
     if (k1 > k2)  
         error(_("Lower band %d > upper band %d"), k1, k2);  
297      for (j = 0; j < n; j++) {      for (j = 0; j < n; j++) {
298          int i1 = j - k2, i2 = j + 1 - k1;          int i1 = j - k2, i2 = j + 1 - k1;
299          for (i = 0; i < i1; i++) ax[i + j * m] = 0.;          for (i = 0; i < i1; i++) ax[i + j * m] = 0.;
# Line 311  Line 315 
315      UNPROTECT(2);      UNPROTECT(2);
316      return aa;      return aa;
317  }  }
318    }
319    
320    SEXP ddense_to_symmetric(SEXP x, SEXP uplo)
321    /* Class of the value will be dsyMatrix
322     */
323    {
324        SEXP xx = PROTECT(dup_mMatrix_as_dgeMatrix(x)),
325            ans = PROTECT(NEW_OBJECT(MAKE_CLASS("dsyMatrix"))), dns;
326    
327        /* Copy xx to ans;
328         * Because slots of xx are freshly allocated and xx will not be
329         * used, we use the slots themselves and don't duplicate */
330        dns = GET_SLOT(xx, Matrix_DimNamesSym);
331        if(!equal_string_vectors(VECTOR_ELT(dns,0),
332                                 VECTOR_ELT(dns,1))) {
333            /* need _symmetric_ dimnames */
334            if(*CHAR(asChar(uplo)) == 'U')
335                SET_VECTOR_ELT(dns,0, VECTOR_ELT(dns,1));
336            else
337                SET_VECTOR_ELT(dns,1, VECTOR_ELT(dns,0));
338        }
339    
340        SET_SLOT(ans, Matrix_xSym,        GET_SLOT(xx, Matrix_xSym));
341        SET_SLOT(ans, Matrix_DimSym,      GET_SLOT(xx, Matrix_DimSym));
342        SET_SLOT(ans, Matrix_DimNamesSym, dns);
343        SET_SLOT(ans, Matrix_uploSym,     ScalarString(asChar(uplo)));
344    
345        UNPROTECT(2);
346        return ans;
347    }
348    
349    SEXP ddense_symmpart(SEXP x)
350    /* Class of the value will be dsyMatrix */
351    {
352        SEXP xx = PROTECT(dup_mMatrix_as_dgeMatrix(x));
353        int *adims = INTEGER(GET_SLOT(xx, Matrix_DimSym)), n = adims[0];
354    
355        if(n != adims[1]) {
356            UNPROTECT(1);
357            error(_("matrix is not square! (symmetric part)"));
358        } else {
359            SEXP ans = PROTECT(NEW_OBJECT(MAKE_CLASS("dsyMatrix"))), dns;
360            double *ax = REAL(GET_SLOT(xx, Matrix_xSym));
361            int i,j;
362    
363            /* only need to assign the *upper* triangle (uplo = "U");
364             * noting that diagonal remains unchanged */
365            for (j = 0; j < n; j++) {
366                for (i = 0; i < j; i++) {
367                    ax[j * n + i] = (ax[j * n + i] + ax[i * n + j]) / 2.;
368                }
369            }
370    
371            /* Copy xx to ans;
372             * Because slots of xx are freshly allocated and xx will not be
373             * used, we use the slots themselves and don't duplicate */
374    
375            dns = GET_SLOT(xx, Matrix_DimNamesSym);
376            if(!equal_string_vectors(VECTOR_ELT(dns,0),
377                                     VECTOR_ELT(dns,1))) {
378                /* need _symmetric_ dimnames */
379                SET_VECTOR_ELT(dns,0, VECTOR_ELT(dns,1));/* ==> uplo = "U" */
380            }
381    
382            SET_SLOT(ans, Matrix_xSym,        GET_SLOT(xx, Matrix_xSym));
383            SET_SLOT(ans, Matrix_DimSym,      GET_SLOT(xx, Matrix_DimSym));
384            SET_SLOT(ans, Matrix_DimNamesSym, dns);
385            SET_SLOT(ans, Matrix_uploSym,     mkString("U"));
386    
387            UNPROTECT(2);
388            return ans;
389        }
390    }
391    
392    SEXP ddense_skewpart(SEXP x)
393    /* Class of the value will be dgeMatrix */
394    {
395        SEXP xx = PROTECT(dup_mMatrix_as_dgeMatrix(x));
396        int *adims = INTEGER(GET_SLOT(xx, Matrix_DimSym)), n = adims[0];
397    
398        if(n != adims[1]) {
399            UNPROTECT(1);
400            error(_("matrix is not square! (skew-symmetric part)"));
401        } else {
402            SEXP ans = PROTECT(NEW_OBJECT(MAKE_CLASS("dgeMatrix"))), dns;
403            double *ax = REAL(GET_SLOT(xx, Matrix_xSym));
404            int i,j;
405    
406            for (j = 0; j < n; j++) {
407                ax[j * n + j] = 0.;
408                for (i = 0; i < j; i++) {
409                    double s = (ax[j * n + i] - ax[i * n + j]) / 2.;
410                    ax[j * n + i] =  s;
411                    ax[i * n + j] = -s;
412                }
413            }
414    
415            /* Copy xx to ans;
416             * Because slots of xx are freshly allocated and xx will not be
417             * used, we use the slots themselves and don't duplicate */
418    
419            dns = GET_SLOT(xx, Matrix_DimNamesSym);
420            if(!equal_string_vectors(VECTOR_ELT(dns,0),
421                                     VECTOR_ELT(dns,1))) {
422                /* need _symmetric_ dimnames */
423                SET_VECTOR_ELT(dns,0, VECTOR_ELT(dns,1));/* uplo = "U" */
424            }
425    
426            SET_SLOT(ans, Matrix_xSym,        GET_SLOT(xx, Matrix_xSym));
427            SET_SLOT(ans, Matrix_DimSym,      GET_SLOT(xx, Matrix_DimSym));
428            SET_SLOT(ans, Matrix_DimNamesSym, dns);
429            SET_SLOT(ans, Matrix_uploSym,     mkString("U"));
430    
431            UNPROTECT(2);
432            return ans;
433        }
434    }

Legend:
Removed from v.2111  
changed lines
  Added in v.2112

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