SCM

SCM Repository

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

Diff of /pkg/src/dtrMatrix.c

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

revision 633, Sun Mar 13 21:01:15 2005 UTC revision 654, Wed Mar 16 16:18:33 2005 UTC
# Line 81  Line 81 
81      return val;      return val;
82  }  }
83    
84  SEXP dtrMatrix_matrix_solve(SEXP a, SEXP b)  SEXP dtrMatrix_matrix_solve(SEXP a, SEXP b, SEXP classed)
85  {  {
86      SEXP val = PROTECT(duplicate(b));      int cl = asLogical(classed);
87      int *Dim = INTEGER(GET_SLOT(a, Matrix_DimSym)),      SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("dgeMatrix")));
88          *bDim = INTEGER(getAttrib(val, R_DimSymbol));      int *adims = INTEGER(GET_SLOT(a, Matrix_DimSym)),
89            *bdims = INTEGER(cl ? GET_SLOT(b, Matrix_DimSym) :
90                             getAttrib(b, R_DimSymbol));
91        int info, n = bdims[0], nrhs = bdims[1];
92        int sz = n * nrhs;
93      double one = 1.0;      double one = 1.0;
94    
95      if (bDim[0] != Dim[1])      if (*adims != *bdims || bdims[1] < 1 || *adims < 1 || *adims != adims[1])
96          error(_("Dimensions of a (%d,%d) and b (%d,%d) do not conform"),          error(_("Dimensions of system to be solved are inconsistent"));
               Dim[0], Dim[1], bDim[0], bDim[1]);  
97      F77_CALL(dtrsm)("L", CHAR(asChar(GET_SLOT(val, Matrix_uploSym))),      F77_CALL(dtrsm)("L", CHAR(asChar(GET_SLOT(val, Matrix_uploSym))),
98                      "N", CHAR(asChar(GET_SLOT(val, Matrix_diagSym))),                      "N", CHAR(asChar(GET_SLOT(val, Matrix_diagSym))),
99                      bDim, bDim+1, &one,                      &n, &nrhs, &one, REAL(GET_SLOT(a, Matrix_xSym)), &n,
100                      REAL(GET_SLOT(a, Matrix_xSym)), Dim,                      Memcpy(REAL(ALLOC_SLOT(val, Matrix_xSym, REALSXP, sz)),
101                      REAL(val), bDim);                             REAL(cl ? GET_SLOT(b, Matrix_xSym):b), sz), &n);
102        UNPROTECT(1);
103        return val;
104    }
105    
106    SEXP dtrMatrix_matrix_mm(SEXP a, SEXP b, SEXP classed, SEXP right)
107    {
108        int cl = asLogical(classed), rt = asLogical(right);
109        SEXP val = PROTECT(NEW_OBJECT(MAKE_CLASS("dgeMatrix")));
110        int *adims = INTEGER(GET_SLOT(a, Matrix_DimSym)),
111            *bdims = INTEGER(cl ? GET_SLOT(b, Matrix_DimSym) :
112                             getAttrib(b, R_DimSymbol)),
113            *cdims = INTEGER(ALLOC_SLOT(val, Matrix_DimSym, INTSXP, 2));
114        int m, n, sz;
115        double one = 1.;
116    
117        if (!cl && !(isReal(b) && isMatrix(b)))
118            error(_("Argument b must be a numeric matrix"));
119        if (adims[0] != adims[1]) error(_("dtrMatrix in \%*\% must be square"));
120        m = rt ? bdims[0] : adims[0];
121        n = rt ? adims[1] : bdims[1];
122        if ((rt && (adims[0] != m)) || (!rt && (bdims[0] != m)))
123                error(_("Matrices are not conformable for multiplication"));
124        if (m < 1 || n < 1)
125            error(_("Matrices with zero extents cannot be multiplied"));
126        cdims[0] = m; cdims[1] = n; sz = m * n;
127        F77_CALL(dtrmm)(rt ? "R" : "L", CHAR(asChar(GET_SLOT(a, Matrix_uploSym))),
128                        "N", CHAR(asChar(GET_SLOT(a, Matrix_diagSym))), &m, &n,
129                        &one, REAL(GET_SLOT(a, Matrix_xSym)), rt ? &n : &m,
130                        Memcpy(REAL(ALLOC_SLOT(val, Matrix_xSym, REALSXP, sz)),
131                               REAL(cl ? GET_SLOT(b, Matrix_xSym) : b), sz),
132                        rt ? &m : &n);
133      UNPROTECT(1);      UNPROTECT(1);
134      return val;      return val;
135  }  }
# Line 148  Line 182 
182      return ret;      return ret;
183  }  }
184    
 SEXP dtrMatrix_dgeMatrix_mm(SEXP a, SEXP b)  
 {  
     int *adims = INTEGER(GET_SLOT(a, Matrix_DimSym)),  
         *bdims = INTEGER(GET_SLOT(b, Matrix_DimSym)),  
         m = adims[0], n = bdims[1], k = adims[1];  
     SEXP val = PROTECT(duplicate(b));  
     double one = 1.;  
   
     if (bdims[0] != k)  
         error(_("Matrices are not conformable for multiplication"));  
     if (m < 1 || n < 1 || k < 1)  
         error(_("Matrices with zero extents cannot be multiplied"));  
     F77_CALL(dtrmm)("L", CHAR(asChar(GET_SLOT(a, Matrix_uploSym))), "N",  
                     CHAR(asChar(GET_SLOT(a, Matrix_diagSym))),  
                     adims, bdims+1, &one,  
                     REAL(GET_SLOT(a, Matrix_xSym)), adims,  
                     REAL(GET_SLOT(val, Matrix_xSym)), bdims);  
     UNPROTECT(1);  
     return val;  
 }  
   
185  SEXP dtrMatrix_dgeMatrix_mm_R(SEXP a, SEXP b)  SEXP dtrMatrix_dgeMatrix_mm_R(SEXP a, SEXP b)
186  {  {
187      int *adims = INTEGER(GET_SLOT(a, Matrix_DimSym)),      int *adims = INTEGER(GET_SLOT(a, Matrix_DimSym)),

Legend:
Removed from v.633  
changed lines
  Added in v.654

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