/[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 971 - (hide annotations)
Wed Feb 14 04:40:49 2007 UTC (12 years, 8 months ago) by ksteube
File MIME type: text/plain
File size: 6661 byte(s)
Had to undo commit to new MPI branch. The changes went into the original and
not the branch. The files committed here are exactly the same as revision 969.


1 jgs 150 /*
2 elspeth 616 ************************************************************
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 jgs 150 */
12 jgs 102
13     /**************************************************************/
14    
15     /* Finley: Mesh */
16    
17     /**************************************************************/
18    
19 jgs 150 /* Author: gross@access.edu.au */
20     /* Version: $Id$ */
21 jgs 102
22     /**************************************************************/
23    
24 jgs 150 #include "Mesh.h"
25 jgs 102 #include "IndexList.h"
26    
27     /**************************************************************/
28    
29     /* returns a reference to the matrix pattern */
30    
31 jgs 150 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 jgs 102 /* 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 jgs 150 if (Finley_noError()) {
49 jgs 102 if (reduce_row_order) {
50     if (reduce_col_order) {
51 jgs 150 out=Paso_SystemMatrixPattern_reference(mesh->ReducedReducedPattern);
52 jgs 102 } else {
53 jgs 150 out=Paso_SystemMatrixPattern_reference(mesh->ReducedFullPattern);
54 jgs 102 }
55     } else {
56     if (reduce_col_order) {
57 jgs 150 out=Paso_SystemMatrixPattern_reference(mesh->FullReducedPattern);
58 jgs 102 } else {
59 jgs 150 out=Paso_SystemMatrixPattern_reference(mesh->FullFullPattern);
60 jgs 102 }
61     }
62     }
63     return out;
64     }
65 jgs 150 Paso_SystemMatrixPattern* Finley_makePattern(Finley_Mesh *mesh,bool_t reduce_row_order, bool_t reduce_col_order) {
66 jgs 102 double time0;
67 jgs 123 dim_t i,n;
68     index_t s,itmp,*rowLabel=NULL,*colLabel=NULL, *ptr=NULL, *index=NULL;
69 jgs 102 Finley_IndexList* index_list=NULL;
70 jgs 150 Finley_resetError();
71 ksteube 971
72 jgs 102 time0=Finley_timer();
73     if (reduce_col_order) {
74     n=mesh->Nodes->reducedNumDegreesOfFreedom;
75     colLabel=mesh->Nodes->reducedDegreeOfFreedom;
76     } else {
77     n=mesh->Nodes->numDegreesOfFreedom;
78     colLabel=mesh->Nodes->degreeOfFreedom;
79     }
80    
81     if (reduce_row_order) {
82     n=mesh->Nodes->reducedNumDegreesOfFreedom;
83     rowLabel=mesh->Nodes->reducedDegreeOfFreedom;
84     } else {
85     n=mesh->Nodes->numDegreesOfFreedom;
86     rowLabel=mesh->Nodes->degreeOfFreedom;
87     }
88    
89     index_list=TMPMEMALLOC(n,Finley_IndexList);
90 jgs 123 ptr=MEMALLOC(n+1,index_t);
91 jgs 102 if (! (Finley_checkPtr(index_list) || Finley_checkPtr(ptr)) ) {
92     #pragma omp parallel private(i,s,itmp)
93     {
94     #pragma omp for schedule(static)
95     for(i=0;i<n;++i) {
96     index_list[i].extension=NULL;
97     index_list[i].n=0;
98     }
99     /* insert contributions from element matrices into colums index index_list: */
100 ksteube 971 Finley_IndexList_insertElements(index_list,mesh->Elements,
101 jgs 102 reduce_row_order,rowLabel,reduce_col_order,colLabel);
102 ksteube 971 Finley_IndexList_insertElements(index_list,mesh->FaceElements,
103 jgs 102 reduce_row_order,rowLabel,reduce_col_order,colLabel);
104 ksteube 971 Finley_IndexList_insertElements(index_list,mesh->ContactElements,
105 jgs 102 reduce_row_order,rowLabel,reduce_col_order,colLabel);
106 ksteube 971 Finley_IndexList_insertElements(index_list,mesh->Points,
107 jgs 102 reduce_row_order,rowLabel,reduce_col_order,colLabel);
108     /* get the number of connections per row */
109     #pragma omp for schedule(static)
110 jgs 123 for(i=0;i<n;++i) {
111     ptr[i]=Finley_IndexList_count(&index_list[i]);
112     }
113 jgs 102 /* accumalate ptr */
114     /* OMP */
115     #pragma omp master
116     {
117     s=0;
118     for(i=0;i<n;++i) {
119     itmp=ptr[i];
120     ptr[i]=s;
121     s+=itmp;
122     }
123     ptr[n]=s;
124     /* fill index */
125 jgs 123 index=MEMALLOC(ptr[n],index_t);
126 jgs 102 }
127     #pragma omp barrier
128     if (! Finley_checkPtr(index)) {
129     #pragma omp for schedule(static)
130     for(i=0;i<n;++i) Finley_IndexList_toArray(&index_list[i],&index[ptr[i]]);
131     }
132     }
133     }
134     /* all done */
135     /* clean up */
136     if (index_list!=NULL) {
137     #pragma omp parallel for private(i)
138     for(i=0;i<n;++i) Finley_IndexList_free(index_list[i].extension);
139     }
140     TMPMEMFREE(index_list);
141 jgs 149 #ifdef Finley_TRACE
142 jgs 102 printf("timing: mesh to matrix pattern: %.4e sec\n",Finley_timer()-time0);
143 jgs 149 #endif
144 jgs 150 if (! Finley_noError()) {
145 jgs 102 MEMFREE(ptr);
146     MEMFREE(index);
147     return NULL;
148     } else {
149 ksteube 971 return Paso_SystemMatrixPattern_alloc(MATRIX_FORMAT_DEFAULT,n,ptr,index);
150 jgs 102 }
151     }
152     /*
153     * $Log$
154 jgs 150 * Revision 1.5 2005/09/15 03:44:22 jgs
155     * Merge of development branch dev-02 back to main trunk on 2005-09-15
156     *
157 jgs 149 * Revision 1.4 2005/09/01 03:31:35 jgs
158     * Merge of development branch dev-02 back to main trunk on 2005-09-01
159     *
160 jgs 150 * Revision 1.3.2.2 2005/09/07 06:26:19 gross
161     * the solver from finley are put into the standalone package paso now
162     *
163 jgs 149 * Revision 1.3.2.1 2005/08/24 02:02:18 gross
164     * timing output switched off. solver output can be swiched through getSolution(verbose=True) now.
165     *
166 jgs 123 * Revision 1.3 2005/07/08 04:07:52 jgs
167     * Merge of development branch back to main trunk on 2005-07-08
168     *
169 jgs 102 * Revision 1.2 2004/12/15 07:08:33 jgs
170     * *** empty log message ***
171     *
172 jgs 123 * Revision 1.1.2.5 2005/06/30 01:53:55 gross
173     * a bug in coloring fixed
174     *
175     * Revision 1.1.2.4 2005/06/29 02:34:51 gross
176     * some changes towards 64 integers in finley
177     *
178 jgs 102 * Revision 1.1.2.3 2004/11/24 01:37:14 gross
179     * some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now
180     *
181     * Revision 1.1.2.1 2004/11/15 01:01:09 gross
182     * and anotther missing file
183     *
184     */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26