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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 432 - (hide annotations)
Fri Jan 13 07:38:54 2006 UTC (13 years, 11 months ago) by gross
File MIME type: text/plain
File size: 3346 byte(s)
some fixes for openmp
1 jgs 150 /* $Id$ */
2    
3     /**************************************************************/
4    
5     /* Paso: SystemMatrixPatternPattern */
6    
7     /**************************************************************/
8    
9     /* Copyrights by ACcESS Australia 2003, 2004, 2005 */
10     /* Author: gross@access.edu.au */
11    
12     /**************************************************************/
13    
14     #include "Paso.h"
15     #include "SystemMatrixPattern.h"
16    
17     /**************************************************************/
18    
19     /* creates SystemMatrixPattern */
20    
21     Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_unrollBlocks(Paso_SystemMatrixPattern* pattern, \
22 gross 415 int type, dim_t row_block_size,dim_t col_block_size) {
23     index_t index_offset_in=(pattern->type & PATTERN_FORMAT_OFFSET1 ? 1:0);
24     index_t index_offset_out=(type & PATTERN_FORMAT_OFFSET1 ? 1:0);
25    
26     if ((pattern->type & PATTERN_FORMAT_SYM) != (type & PATTERN_FORMAT_SYM)) {
27     Paso_setError(TYPE_ERROR,"Paso_SystemMatrixPattern_unrollBlocks: conversion between symmetric and non-symmetric is not implemented yet");
28     return NULL;
29     }
30 jgs 150 Paso_SystemMatrixPattern*out=NULL;
31     index_t *ptr=NULL,*index=NULL,iPtr;
32     dim_t i,j,k,l;
33     Paso_resetError();
34     dim_t block_size=row_block_size*col_block_size;
35     dim_t new_n_ptr=(pattern->n_ptr)*row_block_size;
36     dim_t new_len=(pattern->len)*block_size;
37    
38     ptr=MEMALLOC(new_n_ptr+1,index_t);
39     index=MEMALLOC(new_len,index_t);
40    
41    
42     if (! ( Paso_checkPtr(ptr) || Paso_checkPtr(index) ) ) {
43     #pragma omp parallel
44     {
45     #pragma omp for private(i) schedule(static)
46 gross 415 for (i=0;i<new_n_ptr+1;++i) ptr[i]=index_offset_out;
47 jgs 150
48     #pragma omp master
49 gross 415 ptr[new_n_ptr]=new_len+index_offset_out;
50 jgs 150
51     #pragma omp for private(i,k) schedule(static)
52     for (i=0;i<pattern->n_ptr;++i)
53 gross 415 for (k=0;k<row_block_size;++k) ptr[i*row_block_size+k]=(pattern->ptr[i]-index_offset_in)*block_size+(pattern->ptr[i+1]-pattern->ptr[i])*col_block_size*k+index_offset_out;
54 jgs 150
55     #pragma omp for private(i,iPtr) schedule(static)
56     for (i=0;i<new_n_ptr;++i)
57 gross 432 for (iPtr=ptr[i]-index_offset_out;iPtr<ptr[i+1]-index_offset_out;++iPtr) index[iPtr]=index_offset_out;
58 jgs 150
59     #pragma omp for private(i,j,iPtr,k) schedule(static)
60     for (i=0;i<pattern->n_ptr;++i) {
61 gross 415 for (iPtr=pattern->ptr[i]-index_offset_in;iPtr<pattern->ptr[i+1]-index_offset_in;++iPtr) {
62 jgs 150 for (k=0;k<row_block_size;++k) {
63     for (j=0;j<col_block_size;++j) {
64 gross 425 index[ptr[i*row_block_size+k]-index_offset_out+(iPtr-(pattern->ptr[i]-index_offset_in))*col_block_size+j]=(pattern->index[iPtr]-index_offset_in)*col_block_size+j+index_offset_out;
65 jgs 150 }
66     }
67     }
68     }
69     }
70     /* create return value */
71 gross 415 out=Paso_SystemMatrixPattern_alloc(type,new_n_ptr,ptr,index);
72 jgs 150 }
73     if (! Paso_noError()) {
74     MEMFREE(index);
75     MEMFREE(ptr);
76     }
77     return out;
78     }
79     /*
80     * $Log$
81     * Revision 1.2 2005/09/15 03:44:39 jgs
82     * Merge of development branch dev-02 back to main trunk on 2005-09-15
83     *
84     * Revision 1.1.2.1 2005/09/05 06:29:47 gross
85     * These files have been extracted from finley to define a stand alone libray for iterative
86     * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
87     * has not been tested yet.
88     *
89     *
90     */

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26