/[escript]/trunk/finley/src/Mesh_getPattern.c
ViewVC logotype

Annotation of /trunk/finley/src/Mesh_getPattern.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 102 - (hide annotations)
Wed Dec 15 07:08:39 2004 UTC (14 years, 11 months ago) by jgs
Original Path: trunk/esys2/finley/src/finleyC/Mesh_getPattern.c
File MIME type: text/plain
File size: 5353 byte(s)
*** empty log message ***

1 jgs 102 /* $Id$ */
2    
3     /**************************************************************/
4    
5     /* Finley: Mesh */
6    
7     /**************************************************************/
8    
9     /* Copyrights by ACcESS Australia 2003, 2004 */
10     /* Author: gross@access.edu.au */
11    
12     /**************************************************************/
13    
14     #include "Finley.h"
15     #include "IndexList.h"
16     #include "System.h"
17     #include "Mesh.h"
18    
19     /**************************************************************/
20    
21     /* returns a reference to the matrix pattern */
22    
23     Finley_SystemMatrixPattern* Finley_getPattern(Finley_Mesh *mesh,int reduce_row_order, int reduce_col_order) {
24     Finley_SystemMatrixPattern *out=NULL;
25     Finley_ErrorCode=NO_ERROR;
26     /* make sure that the requested pattern is available */
27     if (reduce_row_order) {
28     if (reduce_col_order) {
29     if (mesh->ReducedReducedPattern==NULL) mesh->ReducedReducedPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
30     } else {
31     if (mesh->ReducedFullPattern==NULL) mesh->ReducedFullPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
32     }
33     } else {
34     if (reduce_col_order) {
35     if (mesh->FullReducedPattern==NULL) mesh->FullReducedPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
36     } else {
37     if (mesh->FullFullPattern==NULL) mesh->FullFullPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
38     }
39     }
40     if (Finley_ErrorCode==NO_ERROR) {
41     if (reduce_row_order) {
42     if (reduce_col_order) {
43     out=Finley_SystemMatrixPattern_reference(mesh->ReducedReducedPattern);
44     } else {
45     out=Finley_SystemMatrixPattern_reference(mesh->ReducedFullPattern);
46     }
47     } else {
48     if (reduce_col_order) {
49     out=Finley_SystemMatrixPattern_reference(mesh->FullReducedPattern);
50     } else {
51     out=Finley_SystemMatrixPattern_reference(mesh->FullFullPattern);
52     }
53     }
54     }
55     return out;
56     }
57     Finley_SystemMatrixPattern* Finley_makePattern(Finley_Mesh *mesh,int reduce_row_order, int reduce_col_order) {
58     double time0;
59     maybelong i,n,s,itmp;
60     maybelong *rowLabel=NULL,*colLabel=NULL, *ptr=NULL, *index=NULL;
61     Finley_IndexList* index_list=NULL;
62     Finley_ErrorCode=NO_ERROR;
63    
64     time0=Finley_timer();
65     if (reduce_col_order) {
66     n=mesh->Nodes->reducedNumDegreesOfFreedom;
67     colLabel=mesh->Nodes->reducedDegreeOfFreedom;
68     } else {
69     n=mesh->Nodes->numDegreesOfFreedom;
70     colLabel=mesh->Nodes->degreeOfFreedom;
71     }
72    
73     if (reduce_row_order) {
74     n=mesh->Nodes->reducedNumDegreesOfFreedom;
75     rowLabel=mesh->Nodes->reducedDegreeOfFreedom;
76     } else {
77     n=mesh->Nodes->numDegreesOfFreedom;
78     rowLabel=mesh->Nodes->degreeOfFreedom;
79     }
80    
81     index_list=TMPMEMALLOC(n,Finley_IndexList);
82     ptr=MEMALLOC(n+1,maybelong);
83     if (! (Finley_checkPtr(index_list) || Finley_checkPtr(ptr)) ) {
84     #pragma omp parallel private(i,s,itmp)
85     {
86     #pragma omp for schedule(static)
87     for(i=0;i<n;++i) {
88     index_list[i].extension=NULL;
89     index_list[i].n=0;
90     }
91     /* insert contributions from element matrices into colums index index_list: */
92     Finley_IndexList_insertElements(index_list,mesh->Elements,
93     reduce_row_order,rowLabel,reduce_col_order,colLabel);
94     Finley_IndexList_insertElements(index_list,mesh->FaceElements,
95     reduce_row_order,rowLabel,reduce_col_order,colLabel);
96     Finley_IndexList_insertElements(index_list,mesh->ContactElements,
97     reduce_row_order,rowLabel,reduce_col_order,colLabel);
98     Finley_IndexList_insertElements(index_list,mesh->Points,
99     reduce_row_order,rowLabel,reduce_col_order,colLabel);
100     /* get the number of connections per row */
101     #pragma omp for schedule(static)
102     for(i=0;i<n;++i) ptr[i]=Finley_IndexList_count(&index_list[i]);
103     /* accumalate ptr */
104     /* OMP */
105     #pragma omp master
106     {
107     s=0;
108     for(i=0;i<n;++i) {
109     itmp=ptr[i];
110     ptr[i]=s;
111     s+=itmp;
112     }
113     ptr[n]=s;
114     /* fill index */
115     index=MEMALLOC(ptr[n],maybelong);
116     }
117     #pragma omp barrier
118     if (! Finley_checkPtr(index)) {
119     #pragma omp for schedule(static)
120     for(i=0;i<n;++i) Finley_IndexList_toArray(&index_list[i],&index[ptr[i]]);
121     }
122     }
123     }
124     /* all done */
125     /* clean up */
126     if (index_list!=NULL) {
127     #pragma omp parallel for private(i)
128     for(i=0;i<n;++i) Finley_IndexList_free(index_list[i].extension);
129     }
130     TMPMEMFREE(index_list);
131     printf("timing: mesh to matrix pattern: %.4e sec\n",Finley_timer()-time0);
132     if (Finley_ErrorCode!=NO_ERROR) {
133     MEMFREE(ptr);
134     MEMFREE(index);
135     return NULL;
136     } else {
137     return Finley_SystemMatrixPattern_alloc(n,ptr,index);
138     }
139     }
140     /*
141     * $Log$
142     * Revision 1.2 2004/12/15 07:08:33 jgs
143     * *** empty log message ***
144     *
145     * Revision 1.1.2.3 2004/11/24 01:37:14 gross
146     * some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now
147     *
148     * Revision 1.1.2.1 2004/11/15 01:01:09 gross
149     * and anotther missing file
150     *
151     */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26