1 /*****************************************************************************
2 * Copyright (c) 2019 FrontISTR Commons
3 * This software is released under the MIT License, see LICENSE.txt
4 *****************************************************************************/
5
6 #include "hecmw_vis_voxel_gen.h"
7
8 #include <stdio.h>
9 #include <stdlib.h>
10 #include "hecmw_vis_mem_util.h"
11 #include "hecmw_vis_comm_util.h"
12 #include "hecmw_vis_ray_trace.h"
13 #include "hecmw_malloc.h"
14
voxel_gen(double range[6],double c_range[2],int nv[3],double * voxel_dxyz,double * voxel_orig_xyz,int * level,int * voxel_n_neighbor_pe,int ** voxel_neighbor_pe,HECMW_Comm VIS_COMM,int vox_on,int display_range_on,double display_range[6])15 void voxel_gen(double range[6], double c_range[2], int nv[3],
16 double *voxel_dxyz, double *voxel_orig_xyz, int *level,
17 int *voxel_n_neighbor_pe, int **voxel_neighbor_pe,
18 HECMW_Comm VIS_COMM, int vox_on, int display_range_on,
19 double display_range[6]) {
20 int my_rank;
21 int pe_size;
22
23 FILE *fp;
24
25 int i, j, k, i1;
26
27 double trange[6], tc_range[2];
28 double *s_range, tmp_range[6];
29 int n_voxel;
30 double dx, dy, dz;
31 int flag, num_neipe, *neipe;
32 double minx, miny, minz, maxx, maxy, maxz;
33 HECMW_Status stat;
34
35 HECMW_Comm_rank(VIS_COMM, &my_rank);
36 HECMW_Comm_size(VIS_COMM, &pe_size);
37
38 s_range = (double *)HECMW_calloc(6 * pe_size, sizeof(double));
39 if (s_range == NULL) HECMW_vis_memory_exit("s_range");
40 if (display_range_on == 0) {
41 HECMW_Allreduce(&range[0], &trange[0], 1, HECMW_DOUBLE, HECMW_MIN,
42 VIS_COMM);
43 HECMW_Allreduce(&range[1], &trange[1], 1, HECMW_DOUBLE, HECMW_MAX,
44 VIS_COMM);
45 HECMW_Allreduce(&range[2], &trange[2], 1, HECMW_DOUBLE, HECMW_MIN,
46 VIS_COMM);
47 HECMW_Allreduce(&range[3], &trange[3], 1, HECMW_DOUBLE, HECMW_MAX,
48 VIS_COMM);
49 HECMW_Allreduce(&range[4], &trange[4], 1, HECMW_DOUBLE, HECMW_MIN,
50 VIS_COMM);
51 HECMW_Allreduce(&range[5], &trange[5], 1, HECMW_DOUBLE, HECMW_MAX,
52 VIS_COMM);
53 } else if (display_range_on == 1) {
54 for (i = 0; i < 6; i++) trange[i] = display_range[i];
55 }
56 if (my_rank != 0) {
57 HECMW_Send(range, 6, HECMW_DOUBLE, 0, 0, VIS_COMM);
58 /* HECMW_Recv(s_range, 6*pe_size, HECMW_DOUBLE, 0, HECMW_ANY_TAG,
59 * VIS_COMM, &stat);
60 */
61 }
62 if (my_rank == 0) {
63 for (i = 0; i < 6; i++) s_range[i] = range[i];
64 for (j = 1; j < pe_size; j++) {
65 HECMW_Recv(tmp_range, 6, HECMW_DOUBLE, j, HECMW_ANY_TAG, VIS_COMM, &stat);
66 for (i = 0; i < 6; i++) s_range[j * 6 + i] = tmp_range[i];
67 }
68 /* for(j=1;j<pe_size;j++)
69 HECMW_Send(s_range, 6*pe_size,HECMW_DOUBLE, j, 0, VIS_COMM);
70 */
71 }
72
73 HECMW_Allreduce(&c_range[0], &tc_range[0], 1, HECMW_DOUBLE, HECMW_MIN,
74 VIS_COMM);
75 HECMW_Allreduce(&c_range[1], &tc_range[1], 1, HECMW_DOUBLE, HECMW_MAX,
76 VIS_COMM);
77 if (my_rank == 0) {
78 if ((fp = fopen("extent.file", "w")) == NULL)
79 HECMW_vis_print_exit("output file name error: extent.file");
80 fprintf(fp, "The range of the whole data field is\n");
81 fprintf(fp, "Minimum x= %lf Maximum x=%lf\n", trange[0], trange[1]);
82 fprintf(fp, "Minimum y= %lf Maximum y=%lf\n", trange[2], trange[3]);
83 fprintf(fp, "Minimum z= %lf Maximum z=%lf\n", trange[4], trange[5]);
84 fprintf(fp, "The range of color component is\n");
85 fprintf(fp, "Minimum color=%lf maximum color=%lf\n", tc_range[0],
86 tc_range[1]);
87 fprintf(stderr, "Minimum color=%lf maximum color=%lf\n", tc_range[0],
88 tc_range[1]);
89 fprintf(fp, "The subrange of each PE\n");
90 for (j = 0; j < pe_size; j++) {
91 fprintf(fp, "PE %d:\n", j);
92 fprintf(fp, "Minimum x= %lf Maximum x=%lf\n", s_range[j * 6 + 0],
93 s_range[j * 6 + 1]);
94 fprintf(fp, "Minimum y= %lf Maximum y=%lf\n", s_range[j * 6 + 2],
95 s_range[j * 6 + 3]);
96 fprintf(fp, "Minimum z= %lf Maximum z=%lf\n", s_range[j * 6 + 4],
97 s_range[j * 6 + 5]);
98 }
99
100 fclose(fp);
101 if (vox_on == 0) {
102 n_voxel = nv[0] * nv[1] * nv[2];
103 /* fprintf(stderr, "the number of voxel is %d\n", n_voxel);
104 */
105 if ((fp = fopen("voxel.file", "w")) == NULL)
106 HECMW_vis_print_exit("output voxel file error: voxel.file");
107 dx = (trange[1] - trange[0]) / nv[0];
108 dy = (trange[3] - trange[2]) / nv[1];
109 dz = (trange[5] - trange[4]) / nv[2];
110 neipe = (int *)HECMW_calloc(pe_size, sizeof(int));
111 if (neipe == NULL) HECMW_vis_memory_exit("in voxel_gen: neipe");
112 if (n_voxel == 0)
113 HECMW_vis_print_exit(
114 "ERROR: HEC-MW-VIS-E1042: n_voxel_x,n_voxel_y, and n_voxel_z "
115 "cannot be less or equal to 0");
116 /* vox->info = (Voxel_info *)HECMW_calloc(n_voxel, sizeof(Voxel_info));
117 if(vox->info==NULL) {
118 fprintf(stderr, "There is no enough memory for vox\n");
119 exit(0);
120 }
121 vox->n_voxel=n_voxel;
122 */
123
124 for (k = 0; k < nv[2]; k++)
125 for (j = 0; j < nv[1]; j++)
126 for (i = 0; i < nv[0]; i++) {
127 fprintf(fp, "%lf %lf %lf %lf %lf %lf\n", trange[0] + i * dx,
128 trange[2] + j * dy, trange[4] + k * dz, dx, dy, dz);
129 voxel_orig_xyz[(k * nv[1] * nv[0] + j * nv[0] + i) * 3] =
130 trange[0] + i * dx;
131 voxel_orig_xyz[(k * nv[1] * nv[0] + j * nv[0] + i) * 3 + 1] =
132 trange[2] + j * dy;
133 voxel_orig_xyz[(k * nv[1] * nv[0] + j * nv[0] + i) * 3 + 2] =
134 trange[4] + k * dz;
135 voxel_dxyz[(k * nv[1] * nv[0] + j * nv[0] + i) * 3] = dx;
136 voxel_dxyz[(k * nv[1] * nv[0] + j * nv[0] + i) * 3 + 1] = dy;
137 voxel_dxyz[(k * nv[1] * nv[0] + j * nv[0] + i) * 3 + 2] = dz;
138
139 /* fprintf(stderr, "i j k=%d %d %d\n", i,j,k);
140 start_i=i-1; end_i=i+1;
141 start_j=j-1; end_j=j+1;
142 start_k=k-1; end_k=k+1;
143 if(start_i<0) start_i=0;
144 if(end_i>nv[0]-1) end_i=nv[0]-1;
145 if(start_j<0) start_j=0;
146 if(end_j>nv[1]-1) end_j=nv[1]-1;
147 if(start_k<0) start_k=0;
148 if(end_k>nv[2]-1) end_k=nv[2]-1;
149 fprintf(stderr, "%d %d %d %d %d %d\n", start_i, end_i, start_j, end_j, start_k,
150 end_k);
151 fprintf(fp, "%d\n", (end_i-start_i+1)*(end_j-start_j+1)*(end_k-start_k+1));
152 for(k1=start_k;k1<=end_k;k1++)
153 for(j1=start_j;j1<=end_j;j1++)
154 for(i1=start_i;i1<=end_i;i1++)
155 fprintf(fp, "%d\n", k1*nv[0]*nv[1]+j1*nv[0]+i1);
156
157 */
158 minx = trange[0] + i * dx;
159 maxx = minx + dx;
160 miny = trange[2] + j * dy;
161 maxy = miny + dy;
162 minz = trange[4] + k * dz;
163 maxz = minz + dz;
164 num_neipe = 0;
165 for (i1 = 0; i1 < pe_size; i1++) {
166 flag = 1;
167 neipe[i1] = -1;
168 if ((s_range[i1 * 6] > maxx + EPSILON * dx) ||
169 (s_range[i1 * 6 + 1] < minx - EPSILON * dx) ||
170 (s_range[i1 * 6 + 2] > maxy + EPSILON * dy) ||
171 (s_range[i1 * 6 + 3] < miny - EPSILON * dy) ||
172 (s_range[i1 * 6 + 4] > maxz + EPSILON * dz) ||
173 (s_range[i1 * 6 + 5] < minz - EPSILON * dz))
174 flag = 0;
175 if (flag == 1) {
176 neipe[num_neipe] = i1;
177 num_neipe++;
178 }
179 }
180 fprintf(fp, "%d\n", num_neipe);
181 voxel_n_neighbor_pe[k * nv[1] * nv[0] + j * nv[0] + i] = num_neipe;
182 if (num_neipe > 0) {
183 for (i1 = 0; i1 < num_neipe; i1++) {
184 fprintf(fp, "%d ", neipe[i1]);
185
186 voxel_neighbor_pe[k * nv[1] * nv[0] + j * nv[0] + i][i1] =
187 neipe[i1];
188
189 if ((i1 % 8) == 7) fprintf(fp, "\n");
190 }
191 fprintf(fp, "\n");
192 }
193 }
194 fclose(fp);
195 /* send vox information to all other PEs */
196 for (i = 1; i < pe_size; i++) {
197 for (j = 0; j < n_voxel; j++) {
198 HECMW_Send(&(voxel_orig_xyz[j * 3]), 1, HECMW_DOUBLE, i, 0, VIS_COMM);
199 HECMW_Send(&(voxel_orig_xyz[j * 3 + 1]), 1, HECMW_DOUBLE, i, 0,
200 VIS_COMM);
201 HECMW_Send(&(voxel_orig_xyz[j * 3 + 2]), 1, HECMW_DOUBLE, i, 0,
202 VIS_COMM);
203 HECMW_Send(&(voxel_dxyz[j * 3]), 1, HECMW_DOUBLE, i, 0, VIS_COMM);
204 HECMW_Send(&(voxel_dxyz[j * 3 + 1]), 1, HECMW_DOUBLE, i, 0, VIS_COMM);
205 HECMW_Send(&(voxel_dxyz[j * 3 + 2]), 1, HECMW_DOUBLE, i, 0, VIS_COMM);
206 HECMW_Send(&(voxel_n_neighbor_pe[j]), 1, HECMW_INT, i, 0, VIS_COMM);
207 if (voxel_n_neighbor_pe[j] > 0)
208 HECMW_Send(voxel_neighbor_pe[j], n_voxel, HECMW_INT, i, 0,
209 VIS_COMM);
210 }
211 }
212 }
213 }
214 if (my_rank != 0) {
215 if (vox_on == 0) {
216 n_voxel = nv[0] * nv[1] * nv[2];
217 if (n_voxel == 0)
218 HECMW_vis_print_exit(
219 "ERROR: HEC-MW-VIS-E1042: n_voxel_x,n_voxel_y, and n_voxel_z "
220 "cannot be less or equal to 0");
221 for (j = 0; j < n_voxel; j++) {
222 HECMW_Recv(&(voxel_orig_xyz[j * 3]), 1, HECMW_DOUBLE, 0, HECMW_ANY_TAG,
223 VIS_COMM, &stat);
224 HECMW_Recv(&(voxel_orig_xyz[j * 3 + 1]), 1, HECMW_DOUBLE, 0,
225 HECMW_ANY_TAG, VIS_COMM, &stat);
226 HECMW_Recv(&(voxel_orig_xyz[j * 3 + 2]), 1, HECMW_DOUBLE, 0,
227 HECMW_ANY_TAG, VIS_COMM, &stat);
228 HECMW_Recv(&(voxel_dxyz[j * 3]), 1, HECMW_DOUBLE, 0, HECMW_ANY_TAG,
229 VIS_COMM, &stat);
230 HECMW_Recv(&(voxel_dxyz[j * 3 + 1]), 1, HECMW_DOUBLE, 0, HECMW_ANY_TAG,
231 VIS_COMM, &stat);
232 HECMW_Recv(&(voxel_dxyz[j * 3 + 2]), 1, HECMW_DOUBLE, 0, HECMW_ANY_TAG,
233 VIS_COMM, &stat);
234 HECMW_Recv(&(voxel_n_neighbor_pe[j]), 1, HECMW_INT, 0, HECMW_ANY_TAG,
235 VIS_COMM, &stat);
236 if (voxel_n_neighbor_pe[j] > 0) {
237 HECMW_Recv(voxel_neighbor_pe[j], n_voxel, HECMW_INT, 0, HECMW_ANY_TAG,
238 VIS_COMM, &stat);
239 }
240 }
241 }
242 }
243
244 return;
245 }
246