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

Contents of /trunk/finley/src/Assemble_getSize.c

Parent Directory Parent Directory | Revision Log Revision Log


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

1 /* $Id$ */
2
3 /**************************************************************/
4
5 /* assemblage routines: */
6
7 /* calculates the minimum distance between two vertices of elements and assigns the value to each */
8 /* quadrature point in element_size */
9
10
11 /**************************************************************/
12
13 /* Copyrights by ACcESS Australia, 2003,2004 */
14 /* author: gross@access.edu.au */
15 /* Version: $Id$ */
16
17 /**************************************************************/
18
19 #include "escript/Data/DataC.h"
20 #include "Finley.h"
21 #include "Assemble.h"
22 #include "NodeFile.h"
23 #include "ElementFile.h"
24 #include "Util.h"
25 #ifdef _OPENMP
26 #include <omp.h>
27 #endif
28
29 /**************************************************************/
30 void Finley_Assemble_getSize(Finley_NodeFile* nodes, Finley_ElementFile* elements, escriptDataC* element_size) {
31
32 double *local_X=NULL,*element_size_array;
33 int e,n0,n1,q,i,node_offset;
34 double d,diff,min_diff;
35
36 if (nodes==NULL || elements==NULL) return;
37 int NVertices=elements->ReferenceElement->Type->numVertices;
38 int NN=elements->ReferenceElement->Type->numNodes;
39 int NS=elements->ReferenceElement->Type->numShapes;
40 int numQuad=elements->ReferenceElement->numQuadNodes;
41 int numDim=nodes->numDim;
42
43 /* set a few more parameters */
44
45 if (getFunctionSpaceType(element_size)==FINLEY_CONTACT_ELEMENTS_2) {
46 node_offset=NN-NS;
47 } else {
48 node_offset=0;
49 }
50
51 /* check the dimensions of element_size */
52
53 if (numDim!=elements->ReferenceElement->Type->numDim) {
54 Finley_ErrorCode=TYPE_ERROR;
55 sprintf(Finley_ErrorMsg,"Gradient: Spatial and element dimension must match.");
56 } else if (! numSamplesEqual(element_size,numQuad,elements->numElements)) {
57 Finley_ErrorCode=TYPE_ERROR;
58 sprintf(Finley_ErrorMsg,"illegal number of samples of element size Data object");
59 } else if (! isDataPointShapeEqual(element_size,0,&(numDim))) {
60 Finley_ErrorCode=TYPE_ERROR;
61 sprintf(Finley_ErrorMsg,"illegal data point shape of element size Data object");
62 } else if (!isExpanded(element_size)) {
63 Finley_ErrorCode=TYPE_ERROR;
64 sprintf(Finley_ErrorMsg,"expanded Data object is expected for element size.");
65 }
66 /* now we can start: */
67
68 if (Finley_ErrorCode==NO_ERROR) {
69 #pragma omp parallel private(local_X)
70 {
71 /* allocation of work arrays */
72 local_X=THREAD_MEMALLOC(NN*numDim,double);
73 if (! Finley_checkPtr(local_X) ) {
74 /* open the element loop */
75 #pragma omp for private(e,min_diff,diff,n0,n1,d,q,i,element_size_array)
76 for(e=0;e<elements->numElements;e++) {
77 /* gather local coordinates of nodes into local_X(numDim,NN): */
78 Finley_Util_Gather_double(NS,&(elements->Nodes[INDEX2(node_offset,e,NN)]),numDim,nodes->Coordinates,local_X);
79 /* calculate minimal differences */
80 min_diff=-1;
81 for (n0=0;n0<NVertices;n0++) {
82 for (n1=n0+1;n1<NVertices;n1++) {
83 diff=0;
84 for (i=0;i<numDim;i++) {
85 d=local_X[INDEX2(i,n0,numDim)]-local_X[INDEX2(i,n1,numDim)];
86 diff+=d*d;
87 }
88 if (min_diff<0) {
89 min_diff=diff;
90 } else {
91 min_diff=MIN(min_diff,diff);
92 }
93 }
94 }
95 min_diff=sqrt(MAX(min_diff,0));
96 /* set all values to min_diff */
97 element_size_array=getSampleData(element_size,e);
98 for (q=0;q<numQuad;q++) element_size_array[q]=min_diff;
99 }
100 }
101 THREAD_MEMFREE(local_X);
102 }
103 }
104 return;
105 }
106 /*
107 * $Log$
108 * Revision 1.4 2004/12/15 07:08:32 jgs
109 * *** empty log message ***
110 *
111 *
112 *
113 */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26