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

Contents of /trunk-mpi-branch/finley/src/Mesh_getPattern.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1120 - (show annotations)
Tue May 1 01:17:46 2007 UTC (13 years, 1 month ago) by ksteube
File MIME type: text/plain
File size: 8691 byte(s)
MPI branch now runs simplesolve.py with both useMPI=yes and useMPI=no

1 /*
2 ************************************************************
3 * Copyright 2006 by ACcESS MNRF *
4 * *
5 * http://www.access.edu.au *
6 * Primary Business: Queensland, Australia *
7 * Licensed under the Open Software License version 3.0 *
8 * http://www.opensource.org/licenses/osl-3.0.php *
9 * *
10 ************************************************************
11 */
12
13 /**************************************************************/
14
15 /* Finley: Mesh */
16
17 /**************************************************************/
18
19 /* Author: gross@access.edu.au */
20 /* Version: $Id$ */
21
22 /**************************************************************/
23
24 #include "Mesh.h"
25 #include "IndexList.h"
26
27 /**************************************************************/
28
29 /* returns a reference to the matrix pattern */
30
31 Paso_SystemMatrixPattern* Finley_getPattern(Finley_Mesh *mesh,bool_t reduce_row_order, bool_t reduce_col_order) {
32 Paso_SystemMatrixPattern *out=NULL;
33 Finley_resetError();
34 /* make sure that the requested pattern is available */
35 if (reduce_row_order) {
36 if (reduce_col_order) {
37 if (mesh->ReducedReducedPattern==NULL) mesh->ReducedReducedPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
38 } else {
39 if (mesh->ReducedFullPattern==NULL) mesh->ReducedFullPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
40 }
41 } else {
42 if (reduce_col_order) {
43 if (mesh->FullReducedPattern==NULL) mesh->FullReducedPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
44 } else {
45 if (mesh->FullFullPattern==NULL) mesh->FullFullPattern=Finley_makePattern(mesh,reduce_row_order,reduce_col_order);
46 }
47 }
48 if (Finley_noError()) {
49 if (reduce_row_order) {
50 if (reduce_col_order) {
51 out=Paso_SystemMatrixPattern_reference(mesh->ReducedReducedPattern);
52 } else {
53 out=Paso_SystemMatrixPattern_reference(mesh->ReducedFullPattern);
54 }
55 } else {
56 if (reduce_col_order) {
57 out=Paso_SystemMatrixPattern_reference(mesh->FullReducedPattern);
58 } else {
59 out=Paso_SystemMatrixPattern_reference(mesh->FullFullPattern);
60 }
61 }
62 }
63 return out;
64 }
65 Paso_SystemMatrixPattern* Finley_makePattern(Finley_Mesh *mesh,bool_t reduce_row_order, bool_t reduce_col_order) {
66 double time0;
67 Paso_SystemMatrixPattern *pattern = NULL;
68 dim_t i,n, numHops;
69 index_t s,itmp,*rowLabel=NULL,*colLabel=NULL, *ptr=NULL, *index=NULL, *hop=NULL;
70 Finley_IndexList* index_list=NULL;
71 Finley_resetError();
72 Finley_NodeDistribution *row_degreeOfFreedomDistribution=NULL;
73 Finley_NodeDistribution *col_degreeOfFreedomDistribution=NULL;
74 Paso_Distribution *row_distribution=NULL;
75 Paso_Distribution *col_distribution=NULL;
76 #ifndef PASO_MPI
77 index_t tmp[2];
78 Paso_MPIInfo *temp_mpi_info = Paso_MPIInfo_alloc( MPI_COMM_WORLD );
79 #endif
80
81 time0=Finley_timer();
82
83 #ifdef PASO_MPI
84 if (reduce_col_order) {
85 n=mesh->Nodes->reducedDegreeOfFreedomDistribution->numLocal;
86 colLabel=mesh->Nodes->reducedDegreeOfFreedom;
87 } else {
88 n=mesh->Nodes->degreeOfFreedomDistribution->numLocal;
89 colLabel=mesh->Nodes->degreeOfFreedom;
90 }
91
92 if (reduce_row_order) {
93 n=mesh->Nodes->reducedDegreeOfFreedomDistribution->numLocal;
94 rowLabel=mesh->Nodes->reducedDegreeOfFreedom;
95 } else {
96 n=mesh->Nodes->degreeOfFreedomDistribution->numLocal;
97 rowLabel=mesh->Nodes->degreeOfFreedom;
98 }
99 #else
100 if (reduce_col_order) {
101 n=mesh->Nodes->reducedNumDegreesOfFreedom;
102 colLabel=mesh->Nodes->reducedDegreeOfFreedom;
103 } else {
104 n=mesh->Nodes->numDegreesOfFreedom;
105 colLabel=mesh->Nodes->degreeOfFreedom;
106 }
107
108 if (reduce_row_order) {
109 n=mesh->Nodes->reducedNumDegreesOfFreedom;
110 rowLabel=mesh->Nodes->reducedDegreeOfFreedom;
111 } else {
112 n=mesh->Nodes->numDegreesOfFreedom;
113 rowLabel=mesh->Nodes->degreeOfFreedom;
114 }
115 #endif
116
117 index_list=TMPMEMALLOC(n,Finley_IndexList);
118 ptr=MEMALLOC(n+1,index_t);
119 if (! (Finley_checkPtr(index_list) || Finley_checkPtr(ptr)) ) {
120
121 #pragma omp parallel private(i,s,itmp)
122 {
123 #pragma omp for schedule(static)
124 for(i=0;i<n;++i) {
125 index_list[i].extension=NULL;
126 index_list[i].n=0;
127 }
128 /* insert contributions from element matrices into colums index index_list: */
129 Finley_IndexList_insertElements(index_list,mesh,mesh->Elements,
130 reduce_row_order,rowLabel,reduce_col_order,colLabel);
131 Finley_IndexList_insertElements(index_list,mesh,mesh->FaceElements,
132 reduce_row_order,rowLabel,reduce_col_order,colLabel);
133 Finley_IndexList_insertElements(index_list,mesh,mesh->ContactElements,
134 reduce_row_order,rowLabel,reduce_col_order,colLabel);
135 Finley_IndexList_insertElements(index_list,mesh,mesh->Points,
136 reduce_row_order,rowLabel,reduce_col_order,colLabel);
137 /* get the number of connections per row */
138 #pragma omp for schedule(static)
139 for(i=0;i<n;++i) {
140 ptr[i]=Finley_IndexList_count(&index_list[i]);
141 }
142 /* accumalate ptr */
143 /* OMP */
144 #pragma omp master
145 {
146 s=0;
147 for(i=0;i<n;++i) {
148 itmp=ptr[i];
149 ptr[i]=s;
150 s+=itmp;
151 }
152 ptr[n]=s;
153 /* fill index */
154 index=MEMALLOC(ptr[n],index_t);
155 }
156 #pragma omp barrier
157 if (! Finley_checkPtr(index)) {
158 #pragma omp for schedule(static)
159 for(i=0;i<n;++i) Finley_IndexList_toArray(&index_list[i],&index[ptr[i]]);
160 }
161 }
162 }
163 /* all done */
164 /* clean up */
165 if (index_list!=NULL) {
166 #pragma omp parallel for private(i)
167 for(i=0;i<n;++i) Finley_IndexList_free(index_list[i].extension);
168 }
169 TMPMEMFREE(index_list);
170 #ifdef Finley_TRACE
171 printf("timing: mesh to matrix pattern: %.4e sec\n",Finley_timer()-time0);
172 #endif
173 if (Finley_noError()) {
174 // Use a getReference method to get the DOF distributions to avoid memory leaks
175 #ifndef PASO_MPI
176 tmp[0] = 0;
177 if (reduce_col_order) {
178 tmp[1] = mesh->Nodes->reducedNumDegreesOfFreedom;
179 } else {
180 tmp[1] = mesh->Nodes->numDegreesOfFreedom;
181 }
182 col_distribution=Paso_Distribution_alloc(temp_mpi_info,tmp,1,0);
183 if (reduce_row_order) {
184 tmp[1] = mesh->Nodes->reducedNumDegreesOfFreedom;
185 } else {
186 tmp[1] = mesh->Nodes->numDegreesOfFreedom;
187 }
188 row_distribution=Paso_Distribution_alloc(temp_mpi_info,tmp,1,0);
189 /* Allocate MPIInfo for serial job */
190 Paso_MPIInfo_dealloc(temp_mpi_info);
191 #else
192 if (reduce_col_order) {
193 col_degreeOfFreedomDistribution = mesh->Nodes->reducedDegreeOfFreedomDistribution;
194 } else {
195 col_degreeOfFreedomDistribution = mesh->Nodes->degreeOfFreedomDistribution;
196 }
197 if (reduce_row_order) {
198 row_degreeOfFreedomDistribution = mesh->Nodes->reducedDegreeOfFreedomDistribution;
199 } else {
200 row_degreeOfFreedomDistribution = mesh->Nodes->degreeOfFreedomDistribution;
201 }
202 row_distribution=Paso_Distribution_alloc(mesh->MPIInfo, row_degreeOfFreedomDistribution->vtxdist,1,0);
203 col_distribution=Paso_Distribution_alloc(mesh->MPIInfo, col_degreeOfFreedomDistribution->vtxdist,1,0);
204 #endif
205
206 if (Finley_noError()) {
207 Paso_SystemMatrixPattern_makeHops(PATTERN_FORMAT_DEFAULT,col_distribution,ptr,index,&numHops,&hop);
208 if (Finley_noError()) {
209 pattern = Paso_SystemMatrixPattern_alloc(PATTERN_FORMAT_DEFAULT,
210 row_distribution,
211 col_distribution,
212 ptr,index,numHops, hop);
213 if (Finley_noError()) {
214 #ifdef PASO_MPI
215 pattern->output_node_distribution=Finley_NodeDistribution_getReference(row_degreeOfFreedomDistribution);
216 pattern->input_node_distribution=Finley_NodeDistribution_getReference(col_degreeOfFreedomDistribution);
217 #endif
218 }
219 }
220 }
221 Paso_Distribution_free(row_distribution);
222 Paso_Distribution_free(col_distribution);
223 }
224 MEMFREE(hop);
225 if (! Finley_noError()) {
226 MEMFREE(ptr);
227 MEMFREE(index);
228 }
229 return pattern;
230 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26