SCM

SCM Repository

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

Diff of /pkg/src/chm_common.c

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

revision 921, Sat Sep 17 14:27:41 2005 UTC revision 922, Sun Sep 18 16:33:54 2005 UTC
# Line 44  Line 44 
44      case 1:      case 1:
45      case 2:      case 2:
46          ans->xtype = CHOLMOD_REAL;          ans->xtype = CHOLMOD_REAL;
47            ans->x = (void *) REAL(GET_SLOT(x, Matrix_xSym));
48          break;          break;
49      case 3:      case 3:
50      case 4:      case 4:
51      case 5:      case 5:
52          ans->xtype = CHOLMOD_PATTERN;          ans->xtype = CHOLMOD_PATTERN;
53            break;
54        case 6:
55        case 7:
56        case 8:
57            ans->xtype = CHOLMOD_COMPLEX;
58            ans->x = (void *) COMPLEX(GET_SLOT(x, Matrix_xSym));
59            break;
60        }
61                                    /* set the stype */
62        switch(ctype % 3) {
63        case 0: ans->stype = 0; break;
64        case 1:
65            ans->stype =
66                (!strcmp(CHAR(asChar(getAttrib(x, Matrix_uploSym))), "U")) ?
67                1 : -1;
68            break;
69        case 2: error("triangular matrices not yet mapped");
70        }
71    
72        return ans;
73    }
74    
75    /**
76     * Copy the contents of a to an appropriate TsparseMatrix object and,
77     * optionally, free a or free both a and its the pointers to its contents.
78     *
79     * @param a matrix to be converted
80     * @param free 0 - don't free a; > 0 cholmod_free a; < 0 Free a
81     *
82     * @return SEXP containing a copy of a
83     */
84    SEXP chm_sparse_to_SEXP(cholmod_sparse *a, int free)
85    {
86        SEXP ans;
87        char *cl;
88        int *dims, nnz = cholmod_nnz(a, &c);
89                                    /* ensure a is sorted and packed */
90        if (!a->sorted || !a->packed) cholmod_sort(a, &c);
91                                    /* determine the class of the result */
92        switch(a->xtype){
93        case CHOLMOD_PATTERN:
94            cl = (a->stype) ? "lsCMatrix" : "lgCMatrix";
95            break;
96        case CHOLMOD_REAL:
97            cl = (a->stype) ? "dsCMatrix" : "dgCMatrix";
98            break;
99        case CHOLMOD_COMPLEX:
100            cl = (a->stype) ? "zsCMatrix" : "zgCMatrix";
101            break;
102        default: error("unknown xtype in cholmod_sparse object");
103        }
104        ans = PROTECT(NEW_OBJECT(MAKE_CLASS(cl)));
105                                    /* allocate and copy common slots */
106        dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
107        dims[0] = a->nrow; dims[1] = a->ncol;
108        Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_pSym, INTSXP, a->ncol + 1)),
109               (int *) a->p, a->ncol + 1);
110        Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_iSym, INTSXP, nnz)),
111               (int *) a->i, nnz);
112                                    /* copy data slot if present */
113        if (a->xtype == CHOLMOD_REAL)
114            Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, nnz)),
115                   (double *) a->x, nnz);
116        if (a->xtype == CHOLMOD_COMPLEX)
117            error("complex sparse matrix code not yet written");
118    /*      Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, nnz)), */
119    /*             (complex *) a->x, nnz); */
120                                    /* set symmetry attributes */
121        if (a->stype)
122            SET_SLOT(ans, Matrix_uploSym,
123                     mkString((a->stype > 0) ? "U" : "L"));
124        if (free > 0) cholmod_free_sparse(&a, &c);
125        if (free < 0) Free(a);
126        UNPROTECT(1);
127        return ans;
128    }
129    
130    /**
131     * Create a cholmod_triplet object with the contents of x.  Note that
132     * the result should *not* be freed with cholmod_triplet_free.  Use
133     * free or Free on the result.
134     *
135     * @param x pointer to an object that inherits from TsparseMatrix
136     *
137     * @return pointer to a cholmod_triplet object that contains pointers
138     * to the slots of x.
139     */
140    cholmod_triplet *as_cholmod_triplet(SEXP x)
141    {
142        cholmod_triplet *ans = (cholmod_triplet*) malloc(sizeof(cholmod_triplet));
143        char *valid[] = {"dgTMatrix", "dsTMatrix", "dtTMatrix",
144                         "lgTMatrix", "lsTMatrix", "ltTMatrix",
145                         "zgTMatrix", "zsTMatrix", "ztTMatrix",
146                         ""};
147        int *dims, ctype = check_class(CHAR(asChar(getAttrib(x, R_ClassSymbol))),
148                                       valid);
149        SEXP islot;
150    
151        if (ctype < 0) error("invalid class of object to as_cholmod_triplet");
152                                    /* characteristics of the system */
153        ans->itype = CHOLMOD_INT;
154        ans->dtype = CHOLMOD_DOUBLE;
155        ans->x = ans->z = (void *) NULL;
156                                    /* dimensions and nzmax */
157        dims = INTEGER(GET_SLOT(x, Matrix_DimSym));
158        ans->nrow = dims[0];
159        ans->ncol = dims[1];
160        islot = GET_SLOT(x, Matrix_iSym);
161        ans->nnz = ans->nzmax = LENGTH(islot);
162                                    /* slots always present */
163        ans->i = (void *) INTEGER(islot);
164        ans->j = (void *) INTEGER(GET_SLOT(x, Matrix_jSym));
165                                    /* set the xtype and any elements */
166        switch(ctype) {
167        case 0:
168        case 1:
169        case 2:
170            ans->xtype = CHOLMOD_REAL;
171          ans->x = (void *) REAL(GET_SLOT(x, Matrix_xSym));          ans->x = (void *) REAL(GET_SLOT(x, Matrix_xSym));
172          break;          break;
173        case 3:
174        case 4:
175        case 5:
176            ans->xtype = CHOLMOD_PATTERN;
177            break;
178      case 6:      case 6:
179      case 7:      case 7:
180      case 8:      case 8:
# Line 72  Line 196 
196      return ans;      return ans;
197  }  }
198    
199    /**
200     * Copy the contents of a to an appropriate TsparseMatrix object and,
201     * optionally, free a or free both a and its the pointers to its contents.
202     *
203     * @param a matrix to be converted
204     * @param free 0 - don't free a; > 0 cholmod_free a; < 0 Free a
205     *
206     * @return SEXP containing a copy of a
207     */
208    SEXP chm_triplet_to_SEXP(cholmod_triplet *a, int free)
209    {
210        SEXP ans;
211        char *cl;
212        int *dims;
213                                    /* determine the class of the result */
214        switch(a->xtype){
215        case CHOLMOD_PATTERN:
216            cl = (a->stype) ? "lsTMatrix" : "lgTMatrix";
217            break;
218        case CHOLMOD_REAL:
219            cl = (a->stype) ? "dsTMatrix" : "dgTMatrix";
220            break;
221        case CHOLMOD_COMPLEX:
222            cl = (a->stype) ? "zsTMatrix" : "zgTMatrix";
223            break;
224        default: error("unknown xtype in cholmod_triplet object");
225        }
226        ans = PROTECT(NEW_OBJECT(MAKE_CLASS(cl)));
227                                    /* allocate and copy common slots */
228        dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
229        dims[0] = a->nrow; dims[1] = a->ncol;
230        Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_iSym, INTSXP, a->nnz)),
231               (int *) a->i, a->nnz);
232        Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_jSym, INTSXP, a->nnz)),
233               (int *) a->j, a->nnz);
234                                    /* copy data slot if present */
235        if (a->xtype == CHOLMOD_REAL)
236            Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, a->nnz)),
237                   (double *) a->x, a->nnz);
238        if (a->xtype == CHOLMOD_COMPLEX)
239            error("complex sparse matrix code not yet written");
240    /*      Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, a->nnz)), */
241    /*             (complex *) a->x, a->nz); */
242                                    /* set symmetry attributes */
243        if (a->stype)
244            SET_SLOT(ans, Matrix_uploSym,
245                     mkString((a->stype > 0) ? "U" : "L"));
246        if (free > 0) cholmod_free_triplet(&a, &c);
247        if (free < 0) Free(a);
248        UNPROTECT(1);
249        return ans;
250    }
251    
252    /**
253     * Create a cholmod_dense object with the contents of x.  Note that
254     * the result should *not* be freed with cholmod_dense_free.  Use
255     * free or Free on the result.
256     *
257     * @param x pointer to an object that inherits from ddenseMatrix
258     *
259     * @return pointer to a cholmod_dense object that contains a pointer
260     * to the contents of x.
261     */
262    cholmod_dense *as_cholmod_dense(SEXP x)
263    {
264        cholmod_dense *ans = (cholmod_dense*) malloc(sizeof(cholmod_dense));
265        char *valid[] = {"matrix", "dgeMatrix", "lgeMatrix", "zgeMatrix",
266                         ""};
267        int *dims, ctype = check_class(CHAR(asChar(getAttrib(x, R_ClassSymbol))),
268                                       valid);
269    
270        if (ctype < 0) error("invalid class of object to as_cholmod_dense");
271                                    /* characteristics of the system */
272        ans->dtype = CHOLMOD_DOUBLE;
273        ans->x = ans->z = (void *) NULL;
274                                    /* dimensions and nzmax */
275        dims = ctype ? INTEGER(GET_SLOT(x, Matrix_DimSym)) :
276            INTEGER(getAttrib(x, R_DimSymbol));
277        ans->d = ans->nrow = dims[0];
278        ans->ncol = dims[1];
279        ans->nzmax = dims[0] * dims[1];
280                                    /* set the xtype and any elements */
281        switch(ctype) {
282        case 0:
283            if (isReal(x)) {
284                ans->xtype = CHOLMOD_REAL;
285                ans->x = (void *) REAL(x);
286            } else if (isLogical(x)){
287                ans->xtype = CHOLMOD_PATTERN;
288                ans->x = (void *) LOGICAL(x);
289            } else if (isComplex(x)) {
290                ans->xtype = CHOLMOD_COMPLEX;
291                ans->x = (void *) COMPLEX(x);
292            } else error("type not supported");
293            break;
294        case 1:
295            ans->xtype = CHOLMOD_REAL;
296            ans->x = (void *) REAL(GET_SLOT(x, Matrix_xSym));
297            break;
298        case 2:
299            ans->xtype = CHOLMOD_PATTERN;
300            ans->x = (void *) LOGICAL(GET_SLOT(x, Matrix_xSym));
301            break;
302        case 3:
303            ans->xtype = CHOLMOD_COMPLEX;
304            ans->x = (void *) COMPLEX(GET_SLOT(x, Matrix_xSym));
305            break;
306        }
307    
308        return ans;
309    }
310    
311    /**
312     * Copy the contents of a to an appropriate denseMatrix object and,
313     * optionally, free a or free both a and its pointer to its contents.
314     *
315     * @param a matrix to be converted
316     * @param free 0 - don't free a; > 0 cholmod_free a; < 0 Free a
317     *
318     * @return SEXP containing a copy of a
319     */
320    SEXP chm_dense_to_SEXP(cholmod_dense *a, int free)
321    {
322        SEXP ans;
323        char *cl;
324        int *dims, ntot;
325                                    /* determine the class of the result */
326        cl = (a->xtype == CHOLMOD_PATTERN) ? "lgeMatrix" :
327            ((a->xtype == CHOLMOD_REAL) ? "dgeMatrix" :
328             ((a->xtype == CHOLMOD_COMPLEX) ? "zgeMatrix" : ""));
329        if (!strlen(cl)) error("unknown xtype");
330    
331        ans = PROTECT(NEW_OBJECT(MAKE_CLASS(cl)));
332                                    /* allocate and copy common slots */
333        dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
334        dims[0] = a->nrow; dims[1] = a->ncol;
335        ntot = dims[0] * dims[1];
336        if (a->d == a->ncol) {      /* copy data slot if present */
337            if (a->xtype == CHOLMOD_REAL)
338                Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, ntot)),
339                       (double *) a->x, ntot);
340            if (a->xtype == CHOLMOD_COMPLEX)
341                error("complex sparse matrix code not yet written");
342    /*      Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, a->nnz)), */
343    /*             (complex *) a->x, a->nz); */
344        } else error("code for cholmod_dense with holes not yet written");
345    
346        if (free > 0) cholmod_free_dense(&a, &c);
347        if (free < 0) Free(a);
348        UNPROTECT(1);
349        return ans;
350    }
351    

Legend:
Removed from v.921  
changed lines
  Added in v.922

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