/[escript]/branches/doubleplusgood/dudley/src/NodeFile_gather.cpp
ViewVC logotype

Contents of /branches/doubleplusgood/dudley/src/NodeFile_gather.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4332 - (show annotations)
Thu Mar 21 04:21:14 2013 UTC (6 years ago) by jfenwick
File size: 7691 byte(s)
like sand though the hourglass
1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2013 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15
16 /************************************************************************************
17 *
18 * Dudley: Mesh: NodeFile
19 * gathers the NodeFile out from the NodeFile in using the entries
20 * in index[0:out->numNodes-1] which are between min_index and max_index (exclusive)
21 * the node index[i]
22 *
23 ************************************************************************************/
24
25 #include "NodeFile.h"
26
27 /************************************************************************************/
28
29 void Dudley_NodeFile_gatherEntries(dim_t n, index_t * index, index_t min_index, index_t max_index,
30 index_t * Id_out, index_t * Id_in,
31 index_t * Tag_out, index_t * Tag_in,
32 index_t * globalDegreesOfFreedom_out, index_t * globalDegreesOfFreedom_in,
33 dim_t numDim, double *Coordinates_out, double *Coordinates_in)
34 {
35 dim_t i;
36 register index_t k;
37 const index_t range = max_index - min_index;
38 const size_t numDim_size = (size_t) numDim * sizeof(double);
39 #pragma omp parallel for private(i,k) schedule(static)
40 for (i = 0; i < n; i++)
41 {
42 k = index[i] - min_index;
43 if ((k >= 0) && (k < range))
44 {
45 Id_out[i] = Id_in[k];
46 Tag_out[i] = Tag_in[k];
47 globalDegreesOfFreedom_out[i] = globalDegreesOfFreedom_in[k];
48 memcpy(&(Coordinates_out[INDEX2(0, i, numDim)]), &(Coordinates_in[INDEX2(0, k, numDim)]), numDim_size);
49 }
50 }
51 }
52
53 void Dudley_NodeFile_gather(index_t * index, Dudley_NodeFile * in, Dudley_NodeFile * out)
54 {
55 index_t min_id, max_id;
56 Dudley_NodeFile_setGlobalIdRange(&min_id, &max_id, in);
57 Dudley_NodeFile_gatherEntries(out->numNodes, index, min_id, max_id,
58 out->Id, in->Id,
59 out->Tag, in->Tag,
60 out->globalDegreesOfFreedom, in->globalDegreesOfFreedom,
61 out->numDim, out->Coordinates, in->Coordinates);
62 }
63
64 void Dudley_NodeFile_gather_global(index_t * index, Dudley_NodeFile * in, Dudley_NodeFile * out)
65 {
66 index_t min_id, max_id, undefined_node;
67 Esys_MPI_rank buffer_rank, *distribution = NULL;
68 index_t *Id_buffer = NULL, *Tag_buffer = NULL, *globalDegreesOfFreedom_buffer = NULL;
69 double *Coordinates_buffer = NULL;
70 dim_t p, buffer_len, n;
71 char error_msg[100];
72 #ifdef ESYS_MPI
73 Esys_MPI_rank dest, source;
74 MPI_Status status;
75 #endif
76
77 /* get the global range of node ids */
78 Dudley_NodeFile_setGlobalIdRange(&min_id, &max_id, in);
79 undefined_node = min_id - 1;
80
81 distribution = new index_t[in->MPIInfo->size + 1];
82
83 if (!Dudley_checkPtr(distribution))
84 {
85 /* distribute the range of node ids */
86 buffer_len = Esys_MPIInfo_setDistribution(in->MPIInfo, min_id, max_id, distribution);
87 /* allocate buffers */
88 Id_buffer = new index_t[buffer_len];
89 Tag_buffer = new index_t[buffer_len];
90 globalDegreesOfFreedom_buffer = new index_t[buffer_len];
91 Coordinates_buffer = new double[buffer_len * out->numDim];
92 if (!(Dudley_checkPtr(Id_buffer) || Dudley_checkPtr(Tag_buffer) ||
93 Dudley_checkPtr(globalDegreesOfFreedom_buffer) || Dudley_checkPtr(Coordinates_buffer)))
94 {
95 /* fill Id_buffer by the undefined_node marker to check if nodes are defined */
96 #pragma omp parallel for private(n) schedule(static)
97 for (n = 0; n < buffer_len; n++)
98 Id_buffer[n] = undefined_node;
99
100 /* fill the buffer by sending portions around in a circle */
101 #ifdef ESYS_MPI
102 dest = Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
103 source = Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
104 #endif
105 buffer_rank = in->MPIInfo->rank;
106 for (p = 0; p < in->MPIInfo->size; ++p)
107 {
108 if (p > 0)
109 { /* the initial send can be skipped */
110 #ifdef ESYS_MPI
111 MPI_Sendrecv_replace(Id_buffer, buffer_len, MPI_INT,
112 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
113 in->MPIInfo->comm, &status);
114 MPI_Sendrecv_replace(Tag_buffer, buffer_len, MPI_INT,
115 dest, in->MPIInfo->msg_tag_counter + 1, source,
116 in->MPIInfo->msg_tag_counter + 1, in->MPIInfo->comm, &status);
117 MPI_Sendrecv_replace(globalDegreesOfFreedom_buffer, buffer_len, MPI_INT, dest,
118 in->MPIInfo->msg_tag_counter + 2, source, in->MPIInfo->msg_tag_counter + 2,
119 in->MPIInfo->comm, &status);
120 MPI_Sendrecv_replace(Coordinates_buffer, buffer_len * out->numDim, MPI_DOUBLE, dest,
121 in->MPIInfo->msg_tag_counter + 3, source, in->MPIInfo->msg_tag_counter + 3,
122 in->MPIInfo->comm, &status);
123 #endif
124 in->MPIInfo->msg_tag_counter += 4;
125 }
126 buffer_rank = Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank - 1);
127 Dudley_NodeFile_scatterEntries(in->numNodes, in->Id,
128 distribution[buffer_rank], distribution[buffer_rank + 1],
129 Id_buffer, in->Id,
130 Tag_buffer, in->Tag,
131 globalDegreesOfFreedom_buffer, in->globalDegreesOfFreedom,
132 out->numDim, Coordinates_buffer, in->Coordinates);
133 }
134 /* now entries are collected from the buffer again by sending the entries around in a circle */
135 #ifdef ESYS_MPI
136 dest = Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank + 1);
137 source = Esys_MPIInfo_mod(in->MPIInfo->size, in->MPIInfo->rank - 1);
138 #endif
139 buffer_rank = in->MPIInfo->rank;
140 for (p = 0; p < in->MPIInfo->size; ++p)
141 {
142 Dudley_NodeFile_gatherEntries(out->numNodes, index,
143 distribution[buffer_rank], distribution[buffer_rank + 1],
144 out->Id, Id_buffer,
145 out->Tag, Tag_buffer,
146 out->globalDegreesOfFreedom, globalDegreesOfFreedom_buffer,
147 out->numDim, out->Coordinates, Coordinates_buffer);
148 if (p < in->MPIInfo->size - 1)
149 { /* the last send can be skipped */
150 #ifdef ESYS_MPI
151 MPI_Sendrecv_replace(Id_buffer, buffer_len, MPI_INT,
152 dest, in->MPIInfo->msg_tag_counter, source, in->MPIInfo->msg_tag_counter,
153 in->MPIInfo->comm, &status);
154 MPI_Sendrecv_replace(Tag_buffer, buffer_len, MPI_INT,
155 dest, in->MPIInfo->msg_tag_counter + 1, source,
156 in->MPIInfo->msg_tag_counter + 1, in->MPIInfo->comm, &status);
157 MPI_Sendrecv_replace(globalDegreesOfFreedom_buffer, buffer_len, MPI_INT, dest,
158 in->MPIInfo->msg_tag_counter + 2, source, in->MPIInfo->msg_tag_counter + 2,
159 in->MPIInfo->comm, &status);
160 MPI_Sendrecv_replace(Coordinates_buffer, buffer_len * out->numDim, MPI_DOUBLE, dest,
161 in->MPIInfo->msg_tag_counter + 3, source, in->MPIInfo->msg_tag_counter + 3,
162 in->MPIInfo->comm, &status);
163 #endif
164 in->MPIInfo->msg_tag_counter += 4;
165 }
166 buffer_rank = Esys_MPIInfo_mod(in->MPIInfo->size, buffer_rank - 1);
167 }
168 /* check if all nodes are set: */
169 #pragma omp parallel for private(n) schedule(static)
170 for (n = 0; n < out->numNodes; ++n)
171 {
172 if (out->Id[n] == undefined_node)
173 {
174 sprintf(error_msg,
175 "Dudley_NodeFile_gather_global: Node id %d at position %d is referenced but is not defined.",
176 out->Id[n], n);
177 Dudley_setError(VALUE_ERROR, error_msg);
178 }
179 }
180
181 }
182 delete[] Id_buffer;
183 delete[] Tag_buffer;
184 delete[] globalDegreesOfFreedom_buffer;
185 delete[] Coordinates_buffer;
186 }
187 delete[] distribution;
188 /* make sure that the error is global */
189 Esys_MPIInfo_noError(in->MPIInfo);
190 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26