/[escript]/trunk/paso/src/UMFPACK.c
ViewVC logotype

Annotation of /trunk/paso/src/UMFPACK.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1149 - (hide annotations)
Thu May 17 00:30:15 2007 UTC (13 years, 11 months ago) by gross
File MIME type: text/plain
File size: 3985 byte(s)
some fixes from the windows version
1 gross 806 /* $Id$ */
2    
3     /*
4     ********************************************************************************
5     * Copyright 2006 by ACcESS MNRF *
6     * *
7     * http://www.access.edu.au *
8     * Primary Business: Queensland, Australia *
9     * Licensed under the Open Software License version 3.0 *
10     * http://www.opensource.org/licenses/osl-3.0.php *
11     ********************************************************************************
12     */
13    
14     /**************************************************************/
15    
16     /* Paso: interface to the Intel UMFPACK library */
17    
18     /**************************************************************/
19    
20     /* Copyrights by ACcESS Australia 2006 */
21     /* Author: gross@@access.edu.au */
22    
23     /**************************************************************/
24    
25     #include "Paso.h"
26     #include "UMFPACK.h"
27     #ifdef _OPENMP
28     #include <omp.h>
29     #endif
30    
31     /**************************************************************/
32    
33     /* free any extra stuff possibly used by the UMFPACK library */
34    
35     void Paso_UMFPACK_free(Paso_SystemMatrix* A) {
36     #ifdef UMFPACK
37     if (A->solver!=NULL) {
38     Paso_UMFPACK_Handler* pt=(Paso_UMFPACK_Handler*)(A->solver);
39     umfpack_di_free_symbolic(&(pt->symbolic));
40     umfpack_di_free_numeric(&(pt->numeric));
41     MEMFREE(pt);
42     A->solver=NULL;
43     }
44     #endif
45     }
46     /* call the solver: */
47    
48     void Paso_UMFPACK(Paso_SystemMatrix* A,
49     double* out,
50     double* in,
51     Paso_Options* options,
52     Paso_Performance* pp) {
53     #ifdef UMFPACK
54     double time0;
55     double control[UMFPACK_CONTROL], info[UMFPACK_INFO];
56     int error = UMFPACK_OK;
57 gross 1149 Paso_UMFPACK_Handler* pt = NULL;
58 gross 806
59     if (! (A->type & (MATRIX_FORMAT_OFFSET1 + MATRIX_FORMAT_BLK1)) ) {
60     Paso_setError(TYPE_ERROR,"Paso_UMFPACK: UMFPACK requires CSR format with index offset 1 and block size 1.");
61     return;
62     }
63     Performance_startMonitor(pp,PERFORMANCE_ALL);
64 gross 1149 pt = (Paso_UMFPACK_Handler *)(A->solver);
65 gross 806 umfpack_di_defaults(control);
66    
67     if (pt==NULL) {
68     int n = A->num_rows;
69     pt=MEMALLOC(1,Paso_UMFPACK_Handler);
70     if (Paso_checkPtr(pt)) return;
71     A->solver=(void*) pt;
72     time0=Paso_timer();
73 gross 1149 /* call LDU symbolic factorization: */
74 gross 806 error=umfpack_di_symbolic(n,n,A->pattern->ptr,A->pattern->index,A->val,&(pt->symbolic),control,info);
75     if (error != UMFPACK_OK) {
76     Paso_setError(VALUE_ERROR,"symbolic factorization failed.");
77     Paso_UMFPACK_free(A);
78     } else {
79 gross 1149 /* call LDU factorization: */
80 gross 806 error= umfpack_di_numeric(A->pattern->ptr,A->pattern->index,A->val,pt->symbolic,&(pt->numeric),control,info);
81     if (error != UMFPACK_OK) {
82     Paso_setError(ZERO_DIVISION_ERROR,"factorization failed. Most likely the matrix is singular.");
83     Paso_UMFPACK_free(A);
84     }
85     if (options->verbose) printf("timing UMFPACK: LDU factorization: %.4e sec.\n",Paso_timer()-time0);
86     }
87     }
88     if (Paso_noError()) {
89     time0=Paso_timer();
90 gross 1149 /* call forward backward substitution: */
91 gross 806 control[UMFPACK_IRSTEP]=2; /* number of refinement steps */
92     error=umfpack_di_solve(UMFPACK_A,A->pattern->ptr,A->pattern->index,A->val,out,in,pt->numeric,control,info);
93     if (options->verbose) printf("timing UMFPACK: solve: %.4e sec\n",Paso_timer()-time0);
94     if (error != UMFPACK_OK) {
95     Paso_setError(VALUE_ERROR,"forward/backward substition failed. Most likely the matrix is singular.");
96     }
97     }
98     Performance_stopMonitor(pp,PERFORMANCE_ALL);
99     #else
100     Paso_setError(SYSTEM_ERROR,"Paso_UMFPACK:UMFPACK is not avialble.");
101     #endif
102    
103     }
104     /*
105     * $Log$
106     *
107     */

  ViewVC Help
Powered by ViewVC 1.1.26