लघु जवाब: बूस्ट के LAPACK
बाइंडिंग का उपयोग न करें, ये घने मैट्रिक्स के लिए डिजाइन किए गए थे, विरल नहीं मेट्रिसेस, UMFPACK
बजाय का उपयोग करें।
लंबे उत्तर: UMFPACK
ए = = बी को हल करने के लिए सबसे अच्छी पुस्तकालयों में से एक है जब ए बड़ा और विचित्र है।
नीचे नमूना कोड (umfpack_simple.c
के आधार पर) है कि एक साधारण A
और b
उत्पन्न करता है और Ax = b
हल करती है।
#include <stdlib.h>
#include <stdio.h>
#include "umfpack.h"
int *Ap;
int *Ai;
double *Ax;
double *b;
double *x;
/* Generates a sparse matrix problem:
A is n x n tridiagonal matrix
A(i,i-1) = -1;
A(i,i) = 3;
A(i,i+1) = -1;
*/
void generate_sparse_matrix_problem(int n){
int i; /* row index */
int nz; /* nonzero index */
int nnz = 2 + 3*(n-2) + 2; /* number of nonzeros*/
int *Ti; /* row indices */
int *Tj; /* col indices */
double *Tx; /* values */
/* Allocate memory for triplet form */
Ti = malloc(sizeof(int)*nnz);
Tj = malloc(sizeof(int)*nnz);
Tx = malloc(sizeof(double)*nnz);
/* Allocate memory for compressed sparse column form */
Ap = malloc(sizeof(int)*(n+1));
Ai = malloc(sizeof(int)*nnz);
Ax = malloc(sizeof(double)*nnz);
/* Allocate memory for rhs and solution vector */
x = malloc(sizeof(double)*n);
b = malloc(sizeof(double)*n);
/* Construct the matrix A*/
nz = 0;
for (i = 0; i < n; i++){
if (i > 0){
Ti[nz] = i;
Tj[nz] = i-1;
Tx[nz] = -1;
nz++;
}
Ti[nz] = i;
Tj[nz] = i;
Tx[nz] = 3;
nz++;
if (i < n-1){
Ti[nz] = i;
Tj[nz] = i+1;
Tx[nz] = -1;
nz++;
}
b[i] = 0;
}
b[0] = 21; b[1] = 1; b[2] = 17;
/* Convert Triplet to Compressed Sparse Column format */
(void) umfpack_di_triplet_to_col(n,n,nnz,Ti,Tj,Tx,Ap,Ai,Ax,NULL);
/* free triplet format */
free(Ti); free(Tj); free(Tx);
}
int main (void)
{
double *null = (double *) NULL ;
int i, n;
void *Symbolic, *Numeric ;
n = 500000;
generate_sparse_matrix_problem(n);
(void) umfpack_di_symbolic (n, n, Ap, Ai, Ax, &Symbolic, null, null);
(void) umfpack_di_numeric (Ap, Ai, Ax, Symbolic, &Numeric, null, null);
umfpack_di_free_symbolic (&Symbolic);
(void) umfpack_di_solve (UMFPACK_A, Ap, Ai, Ax, x, b, Numeric, null, null);
umfpack_di_free_numeric (&Numeric);
for (i = 0 ; i < 10 ; i++) printf ("x [%d] = %g\n", i, x [i]);
free(b); free(x); free(Ax); free(Ai); free(Ap);
return (0);
}
समारोह generate_sparse_matrix_problem
मैट्रिक्स A
और दाएँ हाथ की ओर b
पैदा करता है। मैट्रिक्स पहले तीन गुना रूप में बनाया गया है। वेक्टर टीआई, टीजे, और टीएक्स पूरी तरह से वर्णन करते हैं। ए ट्रिपलेट फॉर्म बनाना आसान है लेकिन कुशल स्पैर मैट्रिक्स विधियों को संपीड़ित स्पैस कॉलम प्रारूप की आवश्यकता होती है। रूपांतरण umfpack_di_triplet_to_col
के साथ किया जाता है।
umfpack_di_symbolic
के साथ एक प्रतीकात्मक कारक किया जाता है। एक स्पैस A
का LU अपघटन umfpack_di_numeric
के साथ किया जाता है। निचले और ऊपरी त्रिकोणीय हल umfpack_di_solve
के साथ किए जाते हैं।
n
के साथ 500,000 के रूप में, मेरी मशीन पर, पूरे कार्यक्रम को चलाने के लिए लगभग एक सेकंड लगता है। वालग्रिंड रिपोर्ट करता है कि 36 9, 3 9, 64 9 बाइट्स (केवल 352 एमबी से थोड़ा अधिक) आवंटित किए गए थे।
नोट करें page ट्रिपलेट (समन्वय) और संपीड़ित प्रारूप में स्पैर मैट्रिस के लिए बूस्ट के समर्थन पर चर्चा करता है। यदि आप चाहें, तो आप इन बूस्ट ऑब्जेक्ट्स को सरल सरणी UMFPACK
में इनपुट के रूप में बदलने के लिए दिनचर्या लिख सकते हैं।
स्पष्ट बाहर इशारा करते हुए एक मैट्रिक्स है कि आकार सरणी 4x10^12 4x10^14 बाइट्स, या एक ही मैट्रिक्स अकेले के लिए 4 से 400 टेराबाइट्स है:
संबंधित सवाल देखें। (जब तक, जैसा कि नीचे बताया गया है, इसका स्पैस) – cyberconte