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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 102 - (hide annotations)
Wed Dec 15 07:08:39 2004 UTC (14 years, 11 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 jgs 82 /* $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 jgs 102 local_X=THREAD_MEMALLOC(NN*numDim,double);
73 jgs 82 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 jgs 102 * Revision 1.4 2004/12/15 07:08:32 jgs
109 jgs 97 * *** empty log message ***
110 jgs 82 *
111 jgs 97 *
112 jgs 82 *
113     */

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26