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 |
/* assemblage routines: integrates data on quadrature points */ |
16 |
|
17 |
/**************************************************************/ |
18 |
|
19 |
/* Copyrights by ACcESS Australia 2003,2004,2005 */ |
20 |
/* Author: gross@access.edu.au */ |
21 |
/* version: $Id$ */ |
22 |
|
23 |
/**************************************************************/ |
24 |
|
25 |
#include "Assemble.h" |
26 |
#include "Util.h" |
27 |
#ifdef _OPENMP |
28 |
#include <omp.h> |
29 |
#endif |
30 |
|
31 |
/**************************************************************/ |
32 |
|
33 |
void Finley_Assemble_integrate(Finley_NodeFile* nodes, Finley_ElementFile* elements,escriptDataC* data,double* out) { |
34 |
type_t datacase; |
35 |
index_t node_offset; |
36 |
bool_t reducedIntegrationOrder; |
37 |
type_t data_type=getFunctionSpaceType(data); |
38 |
dim_t numComps=getDataPointSize(data); |
39 |
Finley_ElementFile_Jacobeans* jac=NULL; |
40 |
Finley_resetError(); |
41 |
if (nodes==NULL || elements==NULL) return; |
42 |
/* set some parameter */ |
43 |
|
44 |
if (data_type==FINLEY_ELEMENTS) { |
45 |
reducedIntegrationOrder=FALSE; |
46 |
} else if (data_type==FINLEY_FACE_ELEMENTS) { |
47 |
reducedIntegrationOrder=FALSE; |
48 |
} else if (data_type==FINLEY_CONTACT_ELEMENTS_1) { |
49 |
reducedIntegrationOrder=FALSE; |
50 |
} else if (data_type==FINLEY_CONTACT_ELEMENTS_2) { |
51 |
reducedIntegrationOrder=FALSE; |
52 |
} else { |
53 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_integrate: integration of data is not possible."); |
54 |
} |
55 |
|
56 |
/* get access to jacobean */ |
57 |
jac=Finley_ElementFile_borrowJacobeans(elements,nodes,FALSE,reducedIntegrationOrder); |
58 |
|
59 |
if (Finley_noError()) { |
60 |
|
61 |
/* check the shape of the data */ |
62 |
if (! numSamplesEqual(data,jac->ReferenceElement->numQuadNodes,elements->numElements)) { |
63 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_integrate: illegal number of samples of integrant kernel Data object"); |
64 |
} |
65 |
|
66 |
|
67 |
/* now we can start */ |
68 |
|
69 |
if (Finley_noError()) { |
70 |
dim_t q,e,i; |
71 |
double *out_local=NULL, rtmp,*data_array=NULL; |
72 |
for (q=0;q<numComps;q++) out[q]=0; |
73 |
#pragma omp parallel private(q,i,rtmp,data_array,out_local) |
74 |
{ |
75 |
out_local=THREAD_MEMALLOC(numComps,double); |
76 |
if (! Finley_checkPtr(out_local) ) { |
77 |
/* initialize local result */ |
78 |
|
79 |
for (i=0;i<numComps;i++) out_local[i]=0; |
80 |
|
81 |
/* open the element loop */ |
82 |
|
83 |
if (isExpanded(data)) { |
84 |
#pragma omp for private(e) schedule(static) |
85 |
for(e=0;e<elements->numElements;e++) { |
86 |
data_array=getSampleData(data,e); |
87 |
for (q=0;q<jac->ReferenceElement->numQuadNodes;q++) { |
88 |
for (i=0;i<numComps;i++) out_local[i]+=data_array[INDEX2(i,q,numComps)]*jac->volume[INDEX2(q,e,jac->ReferenceElement->numQuadNodes)]; |
89 |
} |
90 |
} |
91 |
} else { |
92 |
#pragma omp for private(e) schedule(static) |
93 |
for(e=0;e<elements->numElements;e++) { |
94 |
data_array=getSampleData(data,e); |
95 |
rtmp=0.; |
96 |
for (q=0;q<jac->ReferenceElement->numQuadNodes;q++) rtmp+=jac->volume[INDEX2(q,e,jac->ReferenceElement->numQuadNodes)]; |
97 |
for (i=0;i<numComps;i++) out_local[i]+=data_array[i]*rtmp; |
98 |
} |
99 |
} |
100 |
/* add local results to global result */ |
101 |
#pragma omp critical |
102 |
for (i=0;i<numComps;i++) out[i]+=out_local[i]; |
103 |
} |
104 |
THREAD_MEMFREE(out_local); |
105 |
} |
106 |
} |
107 |
} |
108 |
} |
109 |
|
110 |
/* |
111 |
* $Log$ |
112 |
* Revision 1.6 2005/09/15 03:44:21 jgs |
113 |
* Merge of development branch dev-02 back to main trunk on 2005-09-15 |
114 |
* |
115 |
* Revision 1.5.2.1 2005/09/07 06:26:18 gross |
116 |
* the solver from finley are put into the standalone package paso now |
117 |
* |
118 |
* Revision 1.5 2005/07/08 04:07:48 jgs |
119 |
* Merge of development branch back to main trunk on 2005-07-08 |
120 |
* |
121 |
* Revision 1.4 2004/12/15 07:08:32 jgs |
122 |
* *** empty log message *** |
123 |
* Revision 1.1.1.1.2.2 2005/06/29 02:34:48 gross |
124 |
* some changes towards 64 integers in finley |
125 |
* |
126 |
* Revision 1.1.1.1.2.1 2004/11/24 01:37:12 gross |
127 |
* some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now |
128 |
* |
129 |
* |
130 |
* |
131 |
*/ |