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

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

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

revision 1930 by ksteube, Thu Sep 25 23:11:13 2008 UTC revision 1931 by artak, Mon Oct 27 01:31:28 2008 UTC
# Line 25  Line 25 
25    
26  #include "Paso.h"  #include "Paso.h"
27  #include "UMFPACK.h"  #include "UMFPACK.h"
28    
29  #ifdef _OPENMP  #ifdef _OPENMP
30  #include <omp.h>  #include <omp.h>
31  #endif  #endif
# Line 101  void Paso_UMFPACK(Paso_SystemMatrix* A, Line 102  void Paso_UMFPACK(Paso_SystemMatrix* A,
102  #else  #else
103      Paso_setError(SYSTEM_ERROR,"Paso_UMFPACK:UMFPACK is not avialble.");      Paso_setError(SYSTEM_ERROR,"Paso_UMFPACK:UMFPACK is not avialble.");
104  #endif  #endif
105    
106    }
107    
108    
109    void Paso_UMFPACK1(Paso_SparseMatrix* A,
110                              double* out,
111                              double* in,
112                              bool_t verbose) {
113    #ifdef UMFPACK
114         double time0;
115         double control[UMFPACK_CONTROL], info[UMFPACK_INFO];
116         int error = UMFPACK_OK;
117         Paso_UMFPACK_Handler* pt = NULL;
118        
119         if (! (A->type & (MATRIX_FORMAT_OFFSET1 + MATRIX_FORMAT_BLK1)) ) {
120            Paso_setError(TYPE_ERROR,"Paso_UMFPACK: UMFPACK requires CSR format with index offset 1 and block size 1.");
121            return;
122         }
123         umfpack_di_defaults(control);
124    
125         if (pt==NULL) {
126            int n = A->numRows;
127            pt=MEMALLOC(1,Paso_UMFPACK_Handler);
128            if (Paso_checkPtr(pt)) return;
129            time0=Paso_timer();
130            /* call LDU symbolic factorization: */
131            error=umfpack_di_symbolic(n,n,A->pattern->ptr,A->pattern->index,A->val,&(pt->symbolic),control,info);
132            if (error != UMFPACK_OK) {
133                 Paso_setError(VALUE_ERROR,"symbolic factorization failed.");
134                 MEMFREE(pt);
135            } else {
136                /* call LDU factorization: */
137                error= umfpack_di_numeric(A->pattern->ptr,A->pattern->index,A->val,pt->symbolic,&(pt->numeric),control,info);
138               if (error != UMFPACK_OK) {
139                 Paso_setError(ZERO_DIVISION_ERROR,"factorization failed. Most likely the matrix is singular.");
140                 MEMFREE(pt);
141               }
142               if (verbose) printf("timing UMFPACK: LDU factorization: %.4e sec.\n",Paso_timer()-time0);
143            }
144         }
145         if (Paso_noError())  {
146            time0=Paso_timer();
147            /* call forward backward substitution: */
148            control[UMFPACK_IRSTEP]=2; /* number of refinement steps */
149            error=umfpack_di_solve(UMFPACK_A,A->pattern->ptr,A->pattern->index,A->val,out,in,pt->numeric,control,info);
150            if (verbose) printf("timing UMFPACK: solve: %.4e sec\n",Paso_timer()-time0);
151            if (error != UMFPACK_OK) {
152                  Paso_setError(VALUE_ERROR,"forward/backward substition failed. Most likely the matrix is singular.");
153            }
154         }
155    #else
156        Paso_setError(SYSTEM_ERROR,"Paso_UMFPACK:UMFPACK is not avialble.");
157    #endif
158    
159  }  }

Legend:
Removed from v.1930  
changed lines
  Added in v.1931

  ViewVC Help
Powered by ViewVC 1.1.26