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

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

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

revision 633 by dhawcroft, Thu Mar 23 05:37:00 2006 UTC revision 1811 by ksteube, Thu Sep 25 23:11:13 2008 UTC
# Line 1  Line 1 
 /* $Id$ */  
1    
2  /*  /*******************************************************
3  ********************************************************************************  *
4  *               Copyright   2006 by ACcESS MNRF                                *  * Copyright (c) 2003-2008 by University of Queensland
5  *                                                                              *  * Earth Systems Science Computational Center (ESSCC)
6  *                 http://www.access.edu.au                                     *  * http://www.uq.edu.au/esscc
7  *           Primary Business: Queensland, Australia                            *  *
8  *     Licensed under the Open Software License version 3.0             *  * Primary Business: Queensland, Australia
9  *        http://www.opensource.org/licenses/osl-3.0.php                        *  * 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    
17  /* Paso: SystemMatrixPatternPattern */  /* Paso: SystemMatrixPattern_unrollBlocks */
18    
19  /**************************************************************/  /**************************************************************/
20    
 /* Copyrights by ACcESS Australia 2003, 2004, 2005 */  
21  /* Author: gross@access.edu.au */  /* Author: gross@access.edu.au */
22    
23  /**************************************************************/  /**************************************************************/
24    
 #include "Paso.h"  
25  #include "SystemMatrixPattern.h"  #include "SystemMatrixPattern.h"
26    
27  /**************************************************************/  /**************************************************************/
28    
29  /* creates SystemMatrixPattern  */  /* creates SystemMatrixPattern  */
30    
31  Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_unrollBlocks(Paso_SystemMatrixPattern* pattern, \  Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_unrollBlocks(Paso_SystemMatrixPattern* pattern,
32                                             int type, dim_t row_block_size,dim_t col_block_size) {                                             int type, dim_t output_block_size,dim_t input_block_size) {
   index_t index_offset_in=(pattern->type & PATTERN_FORMAT_OFFSET1 ? 1:0);  
   index_t index_offset_out=(type & PATTERN_FORMAT_OFFSET1 ? 1:0);  
     
   if ((pattern->type & PATTERN_FORMAT_SYM) != (type & PATTERN_FORMAT_SYM)) {  
       Paso_setError(TYPE_ERROR,"Paso_SystemMatrixPattern_unrollBlocks: conversion between symmetric and non-symmetric is not implemented yet");  
       return NULL;  
   }  
33    Paso_SystemMatrixPattern*out=NULL;    Paso_SystemMatrixPattern*out=NULL;
34    index_t *ptr=NULL,*index=NULL,iPtr;    Paso_Pattern *new_mainPattern=NULL, *new_col_couplePattern=NULL, *new_row_couplePattern=NULL;
35    dim_t i,j,k,l;    Paso_Distribution* new_output_distribution=NULL, *new_input_distribution=NULL;
36    Paso_resetError();    Paso_Connector *new_col_connector=NULL, *new_row_connector=NULL;
37    dim_t block_size=row_block_size*col_block_size;  
38    dim_t new_n_ptr=(pattern->n_ptr)*row_block_size;    new_mainPattern=Paso_Pattern_unrollBlocks(pattern->mainPattern,type,output_block_size,input_block_size);
39    dim_t new_len=(pattern->len)*block_size;    new_col_couplePattern=Paso_Pattern_unrollBlocks(pattern->col_couplePattern,type,output_block_size,input_block_size);
40      new_row_couplePattern=Paso_Pattern_unrollBlocks(pattern->row_couplePattern,type,output_block_size,input_block_size);
41    ptr=MEMALLOC(new_n_ptr+1,index_t);    new_output_distribution=Paso_Distribution_alloc(pattern->output_distribution->mpi_info,
42    index=MEMALLOC(new_len,index_t);                                                    pattern->output_distribution->first_component,
43                                                      output_block_size,0);
44      new_input_distribution=Paso_Distribution_alloc(pattern->input_distribution->mpi_info,
45    if (! ( Paso_checkPtr(ptr) || Paso_checkPtr(index) ) )  {                                                    pattern->input_distribution->first_component,
46       #pragma omp parallel                                                    input_block_size,0);
47       {    new_col_connector=Paso_Connector_unroll(pattern->col_connector,input_block_size);
48          #pragma omp for private(i) schedule(static)    new_row_connector=Paso_Connector_unroll(pattern->row_connector,output_block_size);
49          for (i=0;i<new_n_ptr+1;++i) ptr[i]=index_offset_out;    if (Paso_noError()) {
50         out=Paso_SystemMatrixPattern_alloc(type,
51          #pragma omp master                                          new_output_distribution,
52          ptr[new_n_ptr]=new_len+index_offset_out;                                          new_input_distribution,
53                                            new_mainPattern,
54          #pragma omp for private(i,k) schedule(static)                                          new_col_couplePattern,
55          for (i=0;i<pattern->n_ptr;++i)                                          new_row_couplePattern,
56              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;                                          new_col_connector,
57                                                      new_row_connector);
58          #pragma omp for private(i,iPtr) schedule(static)    }
59          for (i=0;i<new_n_ptr;++i)    Paso_Pattern_free(new_mainPattern);
60              for (iPtr=ptr[i]-index_offset_out;iPtr<ptr[i+1]-index_offset_out;++iPtr) index[iPtr]=index_offset_out;    Paso_Pattern_free(new_col_couplePattern);
61      Paso_Pattern_free(new_row_couplePattern);
62          #pragma omp for private(i,j,iPtr,k) schedule(static)    Paso_Distribution_free(new_output_distribution);
63          for (i=0;i<pattern->n_ptr;++i) {    Paso_Distribution_free(new_input_distribution);
64             for (iPtr=pattern->ptr[i]-index_offset_in;iPtr<pattern->ptr[i+1]-index_offset_in;++iPtr)  {    Paso_Connector_free(new_row_connector);
65                for (k=0;k<row_block_size;++k) {    Paso_Connector_free(new_col_connector);
66                   for (j=0;j<col_block_size;++j) {  
67                      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;    if (Paso_noError()) {
68                   }       return out;
69                }    } else {
70             }       Paso_SystemMatrixPattern_free(out);
71          }       return NULL;
      }  
      /* create return value */  
      out=Paso_SystemMatrixPattern_alloc(type,new_n_ptr,ptr,index);  
   }    
   if (! Paso_noError()) {  
      MEMFREE(index);  
      MEMFREE(ptr);  
72    }    }
   return out;  
73  }  }
 /*  
  * $Log$  
  * Revision 1.2  2005/09/15 03:44:39  jgs  
  * Merge of development branch dev-02 back to main trunk on 2005-09-15  
  *  
  * Revision 1.1.2.1  2005/09/05 06:29:47  gross  
  * These files have been extracted from finley to define a stand alone libray for iterative  
  * linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but  
  * has not been tested yet.  
  *  
  *  
  */  

Legend:
Removed from v.633  
changed lines
  Added in v.1811

  ViewVC Help
Powered by ViewVC 1.1.26