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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1149 - (show annotations)
Thu May 17 00:30:15 2007 UTC (12 years, 5 months ago) by gross
File MIME type: text/plain
File size: 3985 byte(s)
some fixes from the windows version
1 /* $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 Paso_UMFPACK_Handler* pt = NULL;
58
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 pt = (Paso_UMFPACK_Handler *)(A->solver);
65 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 /* call LDU symbolic factorization: */
74 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 /* call LDU factorization: */
80 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 /* call forward backward substitution: */
91 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