paso/src/SystemMatrixPattern.c

/* $Id$ */

/**************************************************************/

/* Paso: SystemMatrixPatternPattern */

/**************************************************************/

/* Copyrights by ACcESS Australia 2003, 2004,2005 */
/* Author: gross@access.edu.au */

/**************************************************************/

#include "Paso.h"
#include "SystemMatrixPattern.h"

/**************************************************************/

/* allocates a SystemMatrixPattern  */

Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_alloc(int type, int n_ptr, index_t* ptr,index_t* index) {
  Paso_SystemMatrixPattern*out;
  index_t index_offset=(type & PATTERN_FORMAT_OFFSET1 ? 1:0);
  index_t loc_min_index,loc_max_index,min_index=index_offset,max_index=index_offset-1;
  dim_t i;
  Paso_resetError();

  if (type & PATTERN_FORMAT_SYM) {
    Paso_setError(TYPE_ERROR,"symmetric matrix pattern is not supported yet");
    return NULL;
  }
  #pragma omp parallel private(loc_min_index,loc_max_index,i)
   {
      loc_min_index=index_offset;
      loc_max_index=index_offset-1;
      if (type & PATTERN_FORMAT_OFFSET1) {
         #pragma omp for schedule(static) 
         for (i=0;i<n_ptr;++i) {
             if (ptr[i]<ptr[i+1]) {
               qsort(&(index[ptr[i]-1]),(size_t)(ptr[i+1]-ptr[i]),sizeof(index_t),Paso_comparIndex); 
               loc_min_index=MIN(loc_min_index,index[ptr[i]-1]);
               loc_max_index=MAX(loc_max_index,index[ptr[i+1]-2]);
             }
         }
      } else {
         #pragma omp for schedule(static) 
         for (i=0;i<n_ptr;++i) {
             if (ptr[i]<ptr[i+1]) {
               qsort(&(index[ptr[i]]),(size_t)(ptr[i+1]-ptr[i]),sizeof(index_t),Paso_comparIndex); 
               loc_min_index=MIN(loc_min_index,index[ptr[i]]);
               loc_max_index=MAX(loc_max_index,index[ptr[i+1]-1]);
             }
         }
      }
      #pragma omp critical
      {
         min_index=MIN(loc_min_index,min_index);
         max_index=MAX(loc_max_index,max_index);
      }
  }
  if (min_index<index_offset) {
    Paso_setError(TYPE_ERROR,"Matrix pattern index out of range.");
    return NULL;
  }

  out=MEMALLOC(1,Paso_SystemMatrixPattern);
  if (Paso_checkPtr(out)) return NULL;
  out->n_ptr=n_ptr;
  out->n_index=max_index-index_offset+1;
  out->ptr=ptr;
  out->index=index;
  out->len=out->ptr[out->n_ptr]-index_offset;
  out->reference_counter=1;
  out->type=type;
  #ifdef Paso_TRACE
  printf("Paso_SystemMatrixPattern_dealloc: system matrix pattern as been allocated.\n");
  #endif
  return out;
}

/* returns a reference to in */

Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_reference(Paso_SystemMatrixPattern* in) {
     if (in!=NULL) {
        ++(in->reference_counter);
     }
     return in;
}
  
/* deallocates a SystemMatrixPattern: */

void Paso_SystemMatrixPattern_dealloc(Paso_SystemMatrixPattern* in) {
  if (in!=NULL) {
     in->reference_counter--;
     if (in->reference_counter<=0) {
        MEMFREE(in->ptr);
        MEMFREE(in->index);
        MEMFREE(in);
        #ifdef Paso_TRACE
        printf("Paso_SystemMatrixPattern_dealloc: system matrix pattern as been deallocated.\n");
        #endif
     }
   }
}
/* *************************************************************/

/*  some routines which help to get the matrix pattern from elements: */

/*  this routine is used by qsort called in Paso_SystemMatrixPattern_alloc */

int Paso_comparIndex(const void *index1,const void *index2){
   index_t Iindex1,Iindex2;
   Iindex1=*(index_t*)index1;
   Iindex2=*(index_t*)index2;
   if (Iindex1<Iindex2) {
      return -1;
   } else {
      if (Iindex1>Iindex2) {
         return 1;
      } else {
         return 0;
      }
   }
}
/*
 * $Log$
 * Revision 1.2  2005/09/15 03:44:38  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.
 *
 *
 */
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			/* allocates a SystemMatrixPattern */
20
21	gross	415	Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_alloc(int type, int n_ptr, index_t* ptr,index_t* index) {
22	jgs	150	Paso_SystemMatrixPattern*out;
23	gross	415	index_t index_offset=(type & PATTERN_FORMAT_OFFSET1 ? 1:0);
24			index_t loc_min_index,loc_max_index,min_index=index_offset,max_index=index_offset-1;
25	jgs	150	dim_t i;
26			Paso_resetError();
27
28	gross	415	if (type & PATTERN_FORMAT_SYM) {
29			Paso_setError(TYPE_ERROR,"symmetric matrix pattern is not supported yet");
30			return NULL;
31			}
32	jgs	150	#pragma omp parallel private(loc_min_index,loc_max_index,i)
33	gross	415	{
34			loc_min_index=index_offset;
35			loc_max_index=index_offset-1;
36			if (type & PATTERN_FORMAT_OFFSET1) {
37			#pragma omp for schedule(static)
38			for (i=0;i<n_ptr;++i) {
39			if (ptr[i]<ptr[i+1]) {
40			qsort(&(index[ptr[i]-1]),(size_t)(ptr[i+1]-ptr[i]),sizeof(index_t),Paso_comparIndex);
41			loc_min_index=MIN(loc_min_index,index[ptr[i]-1]);
42	gross	425	loc_max_index=MAX(loc_max_index,index[ptr[i+1]-2]);
43	gross	415	}
44	jgs	150	}
45	gross	415	} else {
46			#pragma omp for schedule(static)
47			for (i=0;i<n_ptr;++i) {
48			if (ptr[i]<ptr[i+1]) {
49			qsort(&(index[ptr[i]]),(size_t)(ptr[i+1]-ptr[i]),sizeof(index_t),Paso_comparIndex);
50			loc_min_index=MIN(loc_min_index,index[ptr[i]]);
51			loc_max_index=MAX(loc_max_index,index[ptr[i+1]-1]);
52			}
53			}
54			}
55			#pragma omp critical
56			{
57			min_index=MIN(loc_min_index,min_index);
58			max_index=MAX(loc_max_index,max_index);
59			}
60	jgs	150	}
61	gross	415	if (min_index<index_offset) {
62	jgs	150	Paso_setError(TYPE_ERROR,"Matrix pattern index out of range.");
63			return NULL;
64			}
65
66			out=MEMALLOC(1,Paso_SystemMatrixPattern);
67			if (Paso_checkPtr(out)) return NULL;
68			out->n_ptr=n_ptr;
69	gross	425	out->n_index=max_index-index_offset+1;
70	jgs	150	out->ptr=ptr;
71			out->index=index;
72	gross	425	out->len=out->ptr[out->n_ptr]-index_offset;
73	jgs	150	out->reference_counter=1;
74	gross	415	out->type=type;
75	jgs	150	#ifdef Paso_TRACE
76			printf("Paso_SystemMatrixPattern_dealloc: system matrix pattern as been allocated.\n");
77			#endif
78			return out;
79			}
80
81			/* returns a reference to in */
82
83			Paso_SystemMatrixPattern* Paso_SystemMatrixPattern_reference(Paso_SystemMatrixPattern* in) {
84			if (in!=NULL) {
85			++(in->reference_counter);
86			}
87			return in;
88			}
89
90			/* deallocates a SystemMatrixPattern: */
91
92			void Paso_SystemMatrixPattern_dealloc(Paso_SystemMatrixPattern* in) {
93			if (in!=NULL) {
94			in->reference_counter--;
95			if (in->reference_counter<=0) {
96			MEMFREE(in->ptr);
97			MEMFREE(in->index);
98			MEMFREE(in);
99			#ifdef Paso_TRACE
100			printf("Paso_SystemMatrixPattern_dealloc: system matrix pattern as been deallocated.\n");
101			#endif
102			}
103			}
104			}
105			/* *************************************************************/
106
107			/* some routines which help to get the matrix pattern from elements: */
108
109			/* this routine is used by qsort called in Paso_SystemMatrixPattern_alloc */
110
111			int Paso_comparIndex(const void index1,const void index2){
112			index_t Iindex1,Iindex2;
113			Iindex1=(index_t)index1;
114			Iindex2=(index_t)index2;
115			if (Iindex1<Iindex2) {
116			return -1;
117			} else {
118			if (Iindex1>Iindex2) {
119			return 1;
120			} else {
121			return 0;
122			}
123			}
124			}
125			/*
126			* $Log$
127			* Revision 1.2 2005/09/15 03:44:38 jgs
128			* Merge of development branch dev-02 back to main trunk on 2005-09-15
129			*
130			* Revision 1.1.2.1 2005/09/05 06:29:47 gross
131			* These files have been extracted from finley to define a stand alone libray for iterative
132			* linear solvers on the ALTIX. main entry through Paso_solve. this version compiles but
133			* has not been tested yet.
134			*
135			*
136			*/
Name	Value
svn:eol-style	native
svn:keywords	Author Date Id Revision