1 |
|
2 |
/******************************************************* |
3 |
* |
4 |
* Copyright (c) 2003-2008 by University of Queensland |
5 |
* Earth Systems Science Computational Center (ESSCC) |
6 |
* http://www.uq.edu.au/esscc |
7 |
* |
8 |
* Primary Business: Queensland, Australia |
9 |
* Licensed under the Open Software License version 3.0 |
10 |
* http://www.opensource.org/licenses/osl-3.0.php |
11 |
* |
12 |
*******************************************************/ |
13 |
|
14 |
|
15 |
/**************************************************************/ |
16 |
|
17 |
/* assemblage routines: calculates the normal vector at quadrature points on face elements */ |
18 |
|
19 |
/**************************************************************/ |
20 |
|
21 |
#include "Assemble.h" |
22 |
#include "Util.h" |
23 |
#ifdef _OPENMP |
24 |
#include <omp.h> |
25 |
#endif |
26 |
|
27 |
/**************************************************************/ |
28 |
|
29 |
|
30 |
void Finley_Assemble_setNormal(Finley_NodeFile* nodes, Finley_ElementFile* elements, escriptDataC* normal) { |
31 |
double *local_X=NULL, *dVdv=NULL,*normal_array; |
32 |
index_t sign,node_offset; |
33 |
Finley_RefElement* reference_element=NULL; |
34 |
dim_t e,q, NN, NS, numDim, numQuad, numDim_local; |
35 |
if (nodes==NULL || elements==NULL) return; |
36 |
NN=elements->ReferenceElement->Type->numNodes; |
37 |
NS=elements->ReferenceElement->Type->numShapes; |
38 |
numDim=nodes->numDim; |
39 |
Finley_resetError(); |
40 |
if (Finley_Assemble_reducedIntegrationOrder(normal)) { |
41 |
reference_element=elements->ReferenceElementReducedOrder; |
42 |
} else { |
43 |
reference_element=elements->ReferenceElement; |
44 |
} |
45 |
numQuad=reference_element->numQuadNodes; |
46 |
numDim_local=reference_element->Type->numDim; |
47 |
|
48 |
/* set some parameter */ |
49 |
|
50 |
if (getFunctionSpaceType(normal)==FINLEY_CONTACT_ELEMENTS_2) { |
51 |
node_offset=NN-NS; |
52 |
sign=-1; |
53 |
} else { |
54 |
node_offset=0; |
55 |
sign=1; |
56 |
} |
57 |
/* check the dimensions of normal */ |
58 |
if (! (numDim==numDim_local || numDim-1==numDim_local)) { |
59 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_setNormal: Cannot calculate normal vector"); |
60 |
} else if (! isDataPointShapeEqual(normal,1,&(numDim))) { |
61 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_setNormal: illegal number of samples of normal Data object"); |
62 |
} else if (! numSamplesEqual(normal,numQuad,elements->numElements)) { |
63 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_setNormal: illegal number of samples of normal Data object"); |
64 |
} else if (! isDataPointShapeEqual(normal,1,&(numDim))) { |
65 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_setNormal: illegal point data shape of normal Data object"); |
66 |
} else if (!isExpanded(normal)) { |
67 |
Finley_setError(TYPE_ERROR,"Finley_Assemble_setNormal: expanded Data object is expected for normal."); |
68 |
} |
69 |
|
70 |
/* now we can start */ |
71 |
if (Finley_noError()) { |
72 |
#pragma omp parallel private(local_X,dVdv) |
73 |
{ |
74 |
local_X=dVdv=NULL; |
75 |
/* allocation of work arrays */ |
76 |
local_X=THREAD_MEMALLOC(NS*numDim,double); |
77 |
dVdv=THREAD_MEMALLOC(numQuad*numDim*numDim_local,double); |
78 |
if (!(Finley_checkPtr(local_X) || Finley_checkPtr(dVdv) ) ) { |
79 |
/* open the element loop */ |
80 |
#pragma omp for private(e,q,normal_array) schedule(static) |
81 |
for(e=0;e<elements->numElements;e++) { |
82 |
/* gather local coordinates of nodes into local_X: */ |
83 |
Finley_Util_Gather_double(NS,&(elements->Nodes[INDEX2(node_offset,e,NN)]),numDim,nodes->Coordinates,local_X); |
84 |
/* calculate dVdv(i,j,q)=local_X(i,n)*DSDv(n,j,q) */ |
85 |
Finley_Util_SmallMatMult(numDim,numDim_local*numQuad,dVdv,NS,local_X,reference_element->dSdv); |
86 |
/* get normalized vector: */ |
87 |
normal_array=getSampleData(normal,e); |
88 |
Finley_NormalVector(numQuad,numDim,numDim_local,dVdv,normal_array); |
89 |
for (q=0;q<numQuad*numDim;q++) normal_array[q]*=sign; |
90 |
} |
91 |
} |
92 |
THREAD_MEMFREE(local_X); |
93 |
THREAD_MEMFREE(dVdv); |
94 |
} |
95 |
} |
96 |
} |
97 |
/* |
98 |
* $Log$ |
99 |
* Revision 1.6 2005/09/15 03:44:21 jgs |
100 |
* Merge of development branch dev-02 back to main trunk on 2005-09-15 |
101 |
* |
102 |
* Revision 1.5.2.1 2005/09/07 06:26:18 gross |
103 |
* the solver from finley are put into the standalone package paso now |
104 |
* |
105 |
* Revision 1.5 2005/07/08 04:07:48 jgs |
106 |
* Merge of development branch back to main trunk on 2005-07-08 |
107 |
* |
108 |
* Revision 1.4 2004/12/15 07:08:32 jgs |
109 |
* *** empty log message *** |
110 |
* Revision 1.1.1.1.2.3 2005/06/29 02:34:48 gross |
111 |
* some changes towards 64 integers in finley |
112 |
* |
113 |
* Revision 1.1.1.1.2.2 2004/11/24 01:37:13 gross |
114 |
* some changes dealing with the integer overflow in memory allocation. Finley solves 4M unknowns now |
115 |
* |
116 |
* |
117 |
* |
118 |
*/ |