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

Parent Directory | Revision Log | Patch

	revision 1386, Thu Aug 17 22:30:40 2006 UTC	revision 1507, Mon Sep 4 18:41:34 2006 UTC
#	Line 2	Line 2
2	#include "dtCMatrix.h"	#include "dtCMatrix.h"
3	#include "cs_utils.h"	#include "cs_utils.h"
4
5	SEXP tsc_validate(SEXP x)	SEXP dtCMatrix_validate(SEXP x)
6	{	{
7	return triangularMatrix_validate(x);	return triangularMatrix_validate(x);
8	/* see ./dsCMatrix.c or ./dtpMatrix.c  on how to do more testing here */	/* see ./dsCMatrix.c or ./dtpMatrix.c  on how to do more testing here */
9	}	}
10
11		#if 0                           /* no longer used */
12	SEXP tsc_to_dgTMatrix(SEXP x)	SEXP tsc_to_dgTMatrix(SEXP x)
13	{	{
14	SEXP ans;	SEXP ans;
#	Line 44	Line 44
44	}	}
45	return ans;	return ans;
46	}	}
47		#endif
48
49	/**	/**
50	* Derive the column pointer vector for the inverse of L from the parent array	* Derive the column pointer vector for the inverse of L from the parent array
#	Line 114	Line 115
115	return ans;	return ans;
116	}	}
117
	#if 0
	SEXP dtCMatrix_solve(SEXP a)
	{
	SEXP ans = PROTECT(NEW_OBJECT(MAKE_CLASS("dtCMatrix")));
	int lo = uplo_P(a)[0] == 'L', unit = diag_P(a)[0] == 'U',
	n = INTEGER(GET_SLOT(a, Matrix_DimSym))[0], nnz,
	*ai = INTEGER(GET_SLOT(a,Matrix_iSym)),
	*ap = INTEGER(GET_SLOT(a, Matrix_pSym)), *bi,
	*bp = INTEGER(ALLOC_SLOT(ans, Matrix_pSym, INTSXP, n + 1));
	int bnz = 10 * ap[n];         /* initial estimate of nnz in b */
	int *ri = Calloc(bnz, int), *ind = Calloc(n, int), j;
	double *ax = REAL(GET_SLOT(a, Matrix_xSym)), *bx;
	SET_SLOT(ans, Matrix_DimSym, duplicate(GET_SLOT(a, Matrix_DimSym)));
	SET_SLOT(ans, Matrix_DimNamesSym,
	duplicate(GET_SLOT(a, Matrix_DimNamesSym)));
	SET_SLOT(ans, Matrix_uploSym, duplicate(GET_SLOT(a, Matrix_uploSym)));
	SET_SLOT(ans, Matrix_diagSym, duplicate(GET_SLOT(a, Matrix_diagSym)));
	if (!(lo && unit))
	error("code for non-unit or upper triangular not yet written");
	/* Initially bp will contain increasing negative values ending at zero. */
	/* Later we add the negative of bp[0] to all values. */
	bp[n] = 0;
	for (j = n - 1; j >= 0; j--) { /* columns in reverse order */
	int i, i1 = ap[j], i2 = ap[j + 1], k, nr;
	if (i1 < i2) AZERO(ind, n);
	for (i = i1; i < i2; i++) {
	ind[ai[i]] = 1;
	for (k = -bp[ai[i] + 1]; k < -bp[ai[i]]; k++) ind[ri[k]] = 1;
	}
	for (k = 0, nr = 0; k < n; k++) if (ind[k]) nr++;
	if ((nr - bp[j + 1]) > bnz) {
	while (nr > (bnz + bp[j + 1])) bnz *= 2;
	ri = Realloc(ri, bnz, int);
	}
	bp[j] = bp[j + 1] - nr;
	for (k = 0, i = -bp[j + 1]; k < n; k++) if (ind[k]) ri[i++] = k;
	}
	bnz = -bp[0];
	bi = INTEGER(ALLOC_SLOT(ans, Matrix_iSym, INTSXP, bnz));
	bx = REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, bnz));
	for (j = 0; j < n; j++) {
	int bpnew = bp[j] + bnz;
	Memcpy(bi + bpnew, ri - bp[j], bp[j + 1] - bp[j]);
	bp[j] = bpnew;
	}
	/* insert code for calculating the actual values here */
	for (j = 0; j < bnz; j++) bx[j] = 1;
	Free(ind); Free(ri);
	UNPROTECT(1);
	return ans;
	}
	#else
118	SEXP dtCMatrix_solve(SEXP a)	SEXP dtCMatrix_solve(SEXP a)
119	{	{
120	SEXP ans = PROTECT(NEW_OBJECT(MAKE_CLASS("dtCMatrix")));	SEXP ans = PROTECT(NEW_OBJECT(MAKE_CLASS("dtCMatrix")));
#	Line 207	Line 153
153	UNPROTECT(1);	UNPROTECT(1);
154	return ans;	return ans;
155	}	}
	#endif
156
157	SEXP dtCMatrix_matrix_solve(SEXP a, SEXP b, SEXP classed)	SEXP dtCMatrix_matrix_solve(SEXP a, SEXP b, SEXP classed)
158	{	{

---

Colored Diff Long Colored Diff Unidiff Context Diff Side by Side

Legend: Removed from v.1386   changed lines   Added in v.1507

---

R-Forge@R-project.org	ViewVC Help
Powered by ViewVC 1.0.0