1 /******************************************************************************
2  * Copyright 1998-2019 Lawrence Livermore National Security, LLC and other
3  * HYPRE Project Developers. See the top-level COPYRIGHT file for details.
4  *
5  * SPDX-License-Identifier: (Apache-2.0 OR MIT)
6  ******************************************************************************/
7 
8 #include "_hypre_parcsr_ls.h"
9 #include "aux_interp.h"
10 
11 /*---------------------------------------------------------------------------
12  * Auxilary routines for the long range interpolation methods.
13  *  Implemented: "standard", "extended", "multipass", "FF"
14  *--------------------------------------------------------------------------*/
15 /* AHB 11/06: Modification of the above original - takes two
16    communication packages and inserts nodes to position expected for
17    OUT_marker
18 
19    offd nodes from comm_pkg take up first chunk of CF_marker_offd, offd
20    nodes from extend_comm_pkg take up the second chunk of CF_marker_offd. */
21 
22 
23 
hypre_alt_insert_new_nodes(hypre_ParCSRCommPkg * comm_pkg,hypre_ParCSRCommPkg * extend_comm_pkg,HYPRE_Int * IN_marker,HYPRE_Int full_off_procNodes,HYPRE_Int * OUT_marker)24 HYPRE_Int hypre_alt_insert_new_nodes(hypre_ParCSRCommPkg *comm_pkg,
25                           hypre_ParCSRCommPkg *extend_comm_pkg,
26                           HYPRE_Int *IN_marker,
27                           HYPRE_Int full_off_procNodes,
28                           HYPRE_Int *OUT_marker)
29 {
30   hypre_ParCSRCommHandle  *comm_handle;
31 
32   HYPRE_Int i, index, shift;
33 
34   HYPRE_Int num_sends, num_recvs;
35 
36   HYPRE_Int *recv_vec_starts;
37 
38   HYPRE_Int e_num_sends;
39 
40   HYPRE_Int *int_buf_data;
41   HYPRE_Int *e_out_marker;
42 
43 
44   num_sends = hypre_ParCSRCommPkgNumSends(comm_pkg);
45   num_recvs =  hypre_ParCSRCommPkgNumRecvs(comm_pkg);
46   recv_vec_starts = hypre_ParCSRCommPkgRecvVecStarts(comm_pkg);
47 
48   e_num_sends = hypre_ParCSRCommPkgNumSends(extend_comm_pkg);
49 
50 
51   index = hypre_max(hypre_ParCSRCommPkgSendMapStart(comm_pkg, num_sends),
52                     hypre_ParCSRCommPkgSendMapStart(extend_comm_pkg, e_num_sends));
53 
54   int_buf_data = hypre_CTAlloc(HYPRE_Int,  index, HYPRE_MEMORY_HOST);
55 
56   /* orig commpkg data*/
57   index = 0;
58 
59   HYPRE_Int begin = hypre_ParCSRCommPkgSendMapStart(comm_pkg, 0);
60   HYPRE_Int end = hypre_ParCSRCommPkgSendMapStart(comm_pkg, num_sends);
61 #ifdef HYPRE_USING_OPENMP
62 #pragma omp parallel for HYPRE_SMP_SCHEDULE
63 #endif
64   for (i = begin; i < end; ++i) {
65      int_buf_data[i - begin] =
66            IN_marker[hypre_ParCSRCommPkgSendMapElmt(comm_pkg, i)];
67   }
68 
69   comm_handle = hypre_ParCSRCommHandleCreate( 11, comm_pkg, int_buf_data,
70                                               OUT_marker);
71 
72   hypre_ParCSRCommHandleDestroy(comm_handle);
73   comm_handle = NULL;
74 
75   /* now do the extend commpkg */
76 
77   /* first we need to shift our position in the OUT_marker */
78   shift = recv_vec_starts[num_recvs];
79   e_out_marker = OUT_marker + shift;
80 
81   index = 0;
82 
83   begin = hypre_ParCSRCommPkgSendMapStart(extend_comm_pkg, 0);
84   end = hypre_ParCSRCommPkgSendMapStart(extend_comm_pkg, e_num_sends);
85 #ifdef HYPRE_USING_OPENMP
86 #pragma omp parallel for HYPRE_SMP_SCHEDULE
87 #endif
88   for (i = begin; i < end; ++i) {
89      int_buf_data[i - begin] =
90            IN_marker[hypre_ParCSRCommPkgSendMapElmt(extend_comm_pkg, i)];
91   }
92 
93   comm_handle = hypre_ParCSRCommHandleCreate( 11, extend_comm_pkg, int_buf_data,
94                                               e_out_marker);
95 
96   hypre_ParCSRCommHandleDestroy(comm_handle);
97   comm_handle = NULL;
98 
99   hypre_TFree(int_buf_data, HYPRE_MEMORY_HOST);
100 
101   return hypre_error_flag;
102 }
103 
hypre_big_insert_new_nodes(hypre_ParCSRCommPkg * comm_pkg,hypre_ParCSRCommPkg * extend_comm_pkg,HYPRE_Int * IN_marker,HYPRE_Int full_off_procNodes,HYPRE_BigInt offset,HYPRE_BigInt * OUT_marker)104 HYPRE_Int hypre_big_insert_new_nodes(hypre_ParCSRCommPkg *comm_pkg,
105                           hypre_ParCSRCommPkg *extend_comm_pkg,
106                           HYPRE_Int *IN_marker,
107                           HYPRE_Int full_off_procNodes,
108                           HYPRE_BigInt offset,
109                           HYPRE_BigInt *OUT_marker)
110 {
111   hypre_ParCSRCommHandle  *comm_handle;
112 
113   HYPRE_Int i, index, shift;
114 
115   HYPRE_Int num_sends, num_recvs;
116 
117   HYPRE_Int *recv_vec_starts;
118 
119   HYPRE_Int e_num_sends;
120 
121   HYPRE_BigInt *int_buf_data;
122   HYPRE_BigInt *e_out_marker;
123 
124 
125   num_sends = hypre_ParCSRCommPkgNumSends(comm_pkg);
126   num_recvs =  hypre_ParCSRCommPkgNumRecvs(comm_pkg);
127   recv_vec_starts = hypre_ParCSRCommPkgRecvVecStarts(comm_pkg);
128 
129   e_num_sends = hypre_ParCSRCommPkgNumSends(extend_comm_pkg);
130 
131 
132   index = hypre_max(hypre_ParCSRCommPkgSendMapStart(comm_pkg, num_sends),
133                     hypre_ParCSRCommPkgSendMapStart(extend_comm_pkg, e_num_sends));
134 
135   int_buf_data = hypre_CTAlloc(HYPRE_BigInt,  index, HYPRE_MEMORY_HOST);
136 
137   /* orig commpkg data*/
138   index = 0;
139 
140   HYPRE_Int begin = hypre_ParCSRCommPkgSendMapStart(comm_pkg, 0);
141   HYPRE_Int end = hypre_ParCSRCommPkgSendMapStart(comm_pkg, num_sends);
142 #ifdef HYPRE_USING_OPENMP
143 #pragma omp parallel for HYPRE_SMP_SCHEDULE
144 #endif
145   for (i = begin; i < end; ++i) {
146      int_buf_data[i - begin] = offset +
147            (HYPRE_BigInt) IN_marker[hypre_ParCSRCommPkgSendMapElmt(comm_pkg, i)];
148   }
149 
150   comm_handle = hypre_ParCSRCommHandleCreate( 21, comm_pkg, int_buf_data,
151                                               OUT_marker);
152 
153   hypre_ParCSRCommHandleDestroy(comm_handle);
154   comm_handle = NULL;
155 
156   /* now do the extend commpkg */
157 
158   /* first we need to shift our position in the OUT_marker */
159   shift = recv_vec_starts[num_recvs];
160   e_out_marker = OUT_marker + shift;
161 
162   index = 0;
163 
164   begin = hypre_ParCSRCommPkgSendMapStart(extend_comm_pkg, 0);
165   end = hypre_ParCSRCommPkgSendMapStart(extend_comm_pkg, e_num_sends);
166 #ifdef HYPRE_USING_OPENMP
167 #pragma omp parallel for HYPRE_SMP_SCHEDULE
168 #endif
169   for (i = begin; i < end; ++i) {
170      int_buf_data[i - begin] = offset +
171            (HYPRE_BigInt) IN_marker[hypre_ParCSRCommPkgSendMapElmt(extend_comm_pkg, i)];
172   }
173 
174   comm_handle = hypre_ParCSRCommHandleCreate( 21, extend_comm_pkg, int_buf_data,
175                                               e_out_marker);
176 
177   hypre_ParCSRCommHandleDestroy(comm_handle);
178   comm_handle = NULL;
179 
180   hypre_TFree(int_buf_data, HYPRE_MEMORY_HOST);
181 
182   return hypre_error_flag;
183 }
184 
185 /* sort for non-ordered arrays */
hypre_ssort(HYPRE_BigInt * data,HYPRE_Int n)186 HYPRE_Int hypre_ssort(HYPRE_BigInt *data, HYPRE_Int n)
187 {
188   HYPRE_Int i,si;
189   HYPRE_Int change = 0;
190 
191   if(n > 0)
192     for(i = n-1; i > 0; i--){
193       si = hypre_index_of_minimum(data,i+1);
194       if(i != si)
195       {
196          hypre_swap_int(data, i, si);
197          change = 1;
198       }
199     }
200   return change;
201 }
202 
203 /* Auxilary function for hypre_ssort */
hypre_index_of_minimum(HYPRE_BigInt * data,HYPRE_Int n)204 HYPRE_Int hypre_index_of_minimum(HYPRE_BigInt *data, HYPRE_Int n)
205 {
206   HYPRE_Int answer;
207   HYPRE_Int i;
208 
209   answer = 0;
210   for(i = 1; i < n; i++)
211     if(data[answer] < data[i])
212       answer = i;
213 
214   return answer;
215 }
216 
hypre_swap_int(HYPRE_BigInt * data,HYPRE_Int a,HYPRE_Int b)217 void hypre_swap_int(HYPRE_BigInt *data, HYPRE_Int a, HYPRE_Int b)
218 {
219   HYPRE_BigInt temp;
220 
221   temp = data[a];
222   data[a] = data[b];
223   data[b] = temp;
224 
225   return;
226 }
227 
228 /* Initialize CF_marker_offd, CF_marker, P_marker, P_marker_offd, tmp */
hypre_initialize_vecs(HYPRE_Int diag_n,HYPRE_Int offd_n,HYPRE_Int * diag_ftc,HYPRE_BigInt * offd_ftc,HYPRE_Int * diag_pm,HYPRE_Int * offd_pm,HYPRE_Int * tmp_CF)229 void hypre_initialize_vecs(HYPRE_Int diag_n, HYPRE_Int offd_n, HYPRE_Int *diag_ftc, HYPRE_BigInt *offd_ftc,
230                            HYPRE_Int *diag_pm, HYPRE_Int *offd_pm, HYPRE_Int *tmp_CF)
231 {
232   HYPRE_Int i;
233 
234   /* Quicker initialization */
235   if(offd_n < diag_n)
236   {
237 #ifdef HYPRE_USING_OPENMP
238 #pragma omp parallel for HYPRE_SMP_SCHEDULE
239 #endif
240     for(i = 0; i < offd_n; i++)
241     {
242       diag_ftc[i] = -1;
243       offd_ftc[i] = -1;
244       tmp_CF[i] = -1;
245       if(diag_pm != NULL)
246       {  diag_pm[i] = -1; }
247       if(offd_pm != NULL)
248       {  offd_pm[i] = -1;}
249     }
250 #ifdef HYPRE_USING_OPENMP
251 #pragma omp parallel for HYPRE_SMP_SCHEDULE
252 #endif
253     for(i = offd_n; i < diag_n; i++)
254     {
255       diag_ftc[i] = -1;
256       if(diag_pm != NULL)
257       {  diag_pm[i] = -1; }
258     }
259   }
260   else
261   {
262 #ifdef HYPRE_USING_OPENMP
263 #pragma omp parallel for HYPRE_SMP_SCHEDULE
264 #endif
265     for(i = 0; i < diag_n; i++)
266     {
267       diag_ftc[i] = -1;
268       offd_ftc[i] = -1;
269       tmp_CF[i] = -1;
270       if(diag_pm != NULL)
271       {  diag_pm[i] = -1;}
272       if(offd_pm != NULL)
273       {  offd_pm[i] = -1;}
274     }
275 #ifdef HYPRE_USING_OPENMP
276 #pragma omp parallel for HYPRE_SMP_SCHEDULE
277 #endif
278     for(i = diag_n; i < offd_n; i++)
279     {
280       offd_ftc[i] = -1;
281       tmp_CF[i] = -1;
282       if(offd_pm != NULL)
283       {  offd_pm[i] = -1;}
284     }
285   }
286   return;
287 }
288 
289 /* Find nodes that are offd and are not contained in original offd
290  * (neighbors of neighbors) */
hypre_new_offd_nodes(HYPRE_BigInt ** found,HYPRE_Int num_cols_A_offd,HYPRE_Int * A_ext_i,HYPRE_BigInt * A_ext_j,HYPRE_Int num_cols_S_offd,HYPRE_BigInt * col_map_offd,HYPRE_BigInt col_1,HYPRE_BigInt col_n,HYPRE_Int * Sop_i,HYPRE_BigInt * Sop_j,HYPRE_Int * CF_marker_offd)291 static HYPRE_Int hypre_new_offd_nodes(HYPRE_BigInt **found, HYPRE_Int num_cols_A_offd,
292        HYPRE_Int *A_ext_i, HYPRE_BigInt *A_ext_j,
293        HYPRE_Int num_cols_S_offd, HYPRE_BigInt *col_map_offd, HYPRE_BigInt col_1,
294        HYPRE_BigInt col_n, HYPRE_Int *Sop_i, HYPRE_BigInt *Sop_j,
295        HYPRE_Int *CF_marker_offd)
296 {
297 #ifdef HYPRE_PROFILE
298   hypre_profile_times[HYPRE_TIMER_ID_RENUMBER_COLIDX] -= hypre_MPI_Wtime();
299 #endif
300 
301   HYPRE_BigInt big_i1, big_k1;
302   HYPRE_Int i, j, kk;
303   HYPRE_Int got_loc, loc_col;
304 
305   /*HYPRE_Int min;*/
306   HYPRE_Int newoff = 0;
307 
308 #ifdef HYPRE_CONCURRENT_HOPSCOTCH
309   hypre_UnorderedBigIntMap col_map_offd_inverse;
310   hypre_UnorderedBigIntMapCreate(&col_map_offd_inverse, 2*num_cols_A_offd, 16*hypre_NumThreads());
311 
312 #pragma omp parallel for HYPRE_SMP_SCHEDULE
313   for (i = 0; i < num_cols_A_offd; i++)
314   {
315      hypre_UnorderedBigIntMapPutIfAbsent(&col_map_offd_inverse, col_map_offd[i], i);
316   }
317 
318   /* Find nodes that will be added to the off diag list */
319   HYPRE_Int size_offP = A_ext_i[num_cols_A_offd];
320   hypre_UnorderedBigIntSet set;
321   hypre_UnorderedBigIntSetCreate(&set, size_offP, 16*hypre_NumThreads());
322 
323 #pragma omp parallel private(i,j,big_i1)
324   {
325 #pragma omp for HYPRE_SMP_SCHEDULE
326     for (i = 0; i < num_cols_A_offd; i++)
327     {
328      if (CF_marker_offd[i] < 0)
329      {
330       for (j = A_ext_i[i]; j < A_ext_i[i+1]; j++)
331       {
332         big_i1 = A_ext_j[j];
333         if(big_i1 < col_1 || big_i1 >= col_n)
334         {
335           if (!hypre_UnorderedBigIntSetContains(&set, big_i1))
336           {
337             HYPRE_Int k = hypre_UnorderedBigIntMapGet(&col_map_offd_inverse, big_i1);
338             if (-1 == k)
339             {
340                hypre_UnorderedBigIntSetPut(&set, big_i1);
341             }
342             else
343             {
344                A_ext_j[j] = -k - 1;
345             }
346           }
347         }
348       }
349       for (j = Sop_i[i]; j < Sop_i[i+1]; j++)
350       {
351         big_i1 = Sop_j[j];
352         if(big_i1 < col_1 || big_i1 >= col_n)
353         {
354           if (!hypre_UnorderedBigIntSetContains(&set, big_i1))
355           {
356             HYPRE_Int k = hypre_UnorderedBigIntMapGet(&col_map_offd_inverse, big_i1);
357             if (-1 == k)
358             {
359                hypre_UnorderedBigIntSetPut(&set, big_i1);
360             }
361             else
362             {
363                Sop_j[j] = -k - 1;
364             }
365           }
366         }
367       }
368      } /* CF_marker_offd[i] < 0 */
369     } /* for each row */
370   } /* omp parallel */
371 
372   hypre_UnorderedBigIntMapDestroy(&col_map_offd_inverse);
373   HYPRE_BigInt *tmp_found = hypre_UnorderedBigIntSetCopyToArray(&set, &newoff);
374   hypre_UnorderedBigIntSetDestroy(&set);
375 
376   /* Put found in monotone increasing order */
377 #ifdef HYPRE_PROFILE
378   hypre_profile_times[HYPRE_TIMER_ID_MERGE] -= hypre_MPI_Wtime();
379 #endif
380 
381   hypre_UnorderedBigIntMap tmp_found_inverse;
382   if (newoff > 0)
383   {
384     hypre_big_sort_and_create_inverse_map(tmp_found, newoff, &tmp_found, &tmp_found_inverse);
385   }
386 
387 #ifdef HYPRE_PROFILE
388   hypre_profile_times[HYPRE_TIMER_ID_MERGE] += hypre_MPI_Wtime();
389 #endif
390 
391   /* Set column indices for Sop and A_ext such that offd nodes are
392    * negatively indexed */
393 #pragma omp parallel for private(kk,big_k1,got_loc,loc_col) HYPRE_SMP_SCHEDULE
394   for(i = 0; i < num_cols_A_offd; i++)
395   {
396    if (CF_marker_offd[i] < 0)
397    {
398      for(kk = Sop_i[i]; kk < Sop_i[i+1]; kk++)
399      {
400        big_k1 = Sop_j[kk];
401        if(big_k1 > -1 && (big_k1 < col_1 || big_k1 >= col_n))
402        {
403          got_loc = hypre_UnorderedBigIntMapGet(&tmp_found_inverse, big_k1);
404          loc_col = got_loc + num_cols_A_offd;
405          Sop_j[kk] = (HYPRE_BigInt)(-loc_col - 1);
406        }
407      }
408      for (kk = A_ext_i[i]; kk < A_ext_i[i+1]; kk++)
409      {
410        big_k1 = A_ext_j[kk];
411        if(big_k1 > -1 && (big_k1 < col_1 || big_k1 >= col_n))
412        {
413          got_loc = hypre_UnorderedBigIntMapGet(&tmp_found_inverse, big_k1);
414          loc_col = got_loc + num_cols_A_offd;
415          A_ext_j[kk] = (HYPRE_BigInt)(-loc_col - 1);
416        }
417      }
418    }
419   }
420   if (newoff)
421   {
422     hypre_UnorderedBigIntMapDestroy(&tmp_found_inverse);
423   }
424 #else /* !HYPRE_CONCURRENT_HOPSCOTCH */
425   HYPRE_Int size_offP;
426 
427   HYPRE_BigInt *tmp_found;
428   HYPRE_Int min;
429   HYPRE_Int ifound;
430 
431   size_offP = A_ext_i[num_cols_A_offd]+Sop_i[num_cols_A_offd];
432   tmp_found = hypre_CTAlloc(HYPRE_BigInt, size_offP, HYPRE_MEMORY_HOST);
433 
434   /* Find nodes that will be added to the off diag list */
435   for (i = 0; i < num_cols_A_offd; i++)
436   {
437    if (CF_marker_offd[i] < 0)
438    {
439     for (j = A_ext_i[i]; j < A_ext_i[i+1]; j++)
440     {
441       big_i1 = A_ext_j[j];
442       if(big_i1 < col_1 || big_i1 >= col_n)
443       {
444          ifound = hypre_BigBinarySearch(col_map_offd,big_i1,num_cols_A_offd);
445          if(ifound == -1)
446          {
447             tmp_found[newoff]=big_i1;
448             newoff++;
449          }
450          else
451          {
452             A_ext_j[j] = (HYPRE_BigInt)(-ifound-1);
453          }
454       }
455     }
456     for (j = Sop_i[i]; j < Sop_i[i+1]; j++)
457     {
458       big_i1 = Sop_j[j];
459       if(big_i1 < col_1 || big_i1 >= col_n)
460       {
461          ifound = hypre_BigBinarySearch(col_map_offd,big_i1,num_cols_A_offd);
462          if(ifound == -1)
463          {
464             tmp_found[newoff]=big_i1;
465             newoff++;
466          }
467          else
468          {
469             Sop_j[j] = (HYPRE_BigInt)(-ifound-1);
470          }
471       }
472     }
473    }
474   }
475   /* Put found in monotone increasing order */
476   if (newoff > 0)
477   {
478      hypre_BigQsort0(tmp_found,0,newoff-1);
479      ifound = tmp_found[0];
480      min = 1;
481      for (i=1; i < newoff; i++)
482      {
483        if (tmp_found[i] > ifound)
484        {
485           ifound = tmp_found[i];
486           tmp_found[min++] = ifound;
487        }
488      }
489      newoff = min;
490   }
491 
492   /* Set column indices for Sop and A_ext such that offd nodes are
493    * negatively indexed */
494   for(i = 0; i < num_cols_A_offd; i++)
495   {
496    if (CF_marker_offd[i] < 0)
497    {
498      for(kk = Sop_i[i]; kk < Sop_i[i+1]; kk++)
499      {
500        big_k1 = Sop_j[kk];
501        if(big_k1 > -1 && (big_k1 < col_1 || big_k1 >= col_n))
502        {
503           got_loc = hypre_BigBinarySearch(tmp_found,big_k1,newoff);
504           if(got_loc > -1)
505              loc_col = got_loc + num_cols_A_offd;
506           Sop_j[kk] = (HYPRE_BigInt)(-loc_col - 1);
507        }
508      }
509      for (kk = A_ext_i[i]; kk < A_ext_i[i+1]; kk++)
510      {
511        big_k1 = A_ext_j[kk];
512        if(big_k1 > -1 && (big_k1 < col_1 || big_k1 >= col_n))
513        {
514           got_loc = hypre_BigBinarySearch(tmp_found,big_k1,newoff);
515           loc_col = got_loc + num_cols_A_offd;
516           A_ext_j[kk] = (HYPRE_BigInt)(-loc_col - 1);
517        }
518      }
519    }
520   }
521 #endif /* !HYPRE_CONCURRENT_HOPSCOTCH */
522 
523   *found = tmp_found;
524 
525 #ifdef HYPRE_PROFILE
526   hypre_profile_times[HYPRE_TIMER_ID_RENUMBER_COLIDX] += hypre_MPI_Wtime();
527 #endif
528 
529   return newoff;
530 }
531 
hypre_exchange_marker(hypre_ParCSRCommPkg * comm_pkg,HYPRE_Int * IN_marker,HYPRE_Int * OUT_marker)532 HYPRE_Int hypre_exchange_marker(hypre_ParCSRCommPkg *comm_pkg,
533                           HYPRE_Int *IN_marker,
534                           HYPRE_Int *OUT_marker)
535 {
536   HYPRE_Int num_sends = hypre_ParCSRCommPkgNumSends(comm_pkg);
537   HYPRE_Int begin = hypre_ParCSRCommPkgSendMapStart(comm_pkg, 0);
538   HYPRE_Int end = hypre_ParCSRCommPkgSendMapStart(comm_pkg, num_sends);
539   HYPRE_Int *int_buf_data = hypre_CTAlloc(HYPRE_Int, end, HYPRE_MEMORY_HOST);
540 
541   HYPRE_Int i;
542 #ifdef HYPRE_USING_OPENMP
543 #pragma omp parallel for HYPRE_SMP_SCHEDULE
544 #endif
545   for (i = begin; i < end; ++i) {
546      int_buf_data[i - begin] =
547            IN_marker[hypre_ParCSRCommPkgSendMapElmt(comm_pkg, i)];
548   }
549 
550   hypre_ParCSRCommHandle *comm_handle = hypre_ParCSRCommHandleCreate( 11, comm_pkg, int_buf_data,
551                                                                       OUT_marker);
552 
553   hypre_ParCSRCommHandleDestroy(comm_handle);
554   hypre_TFree(int_buf_data, HYPRE_MEMORY_HOST);
555 
556   return hypre_error_flag;
557 }
558 
hypre_exchange_interp_data(HYPRE_Int ** CF_marker_offd,HYPRE_Int ** dof_func_offd,hypre_CSRMatrix ** A_ext,HYPRE_Int * full_off_procNodes,hypre_CSRMatrix ** Sop,hypre_ParCSRCommPkg ** extend_comm_pkg,hypre_ParCSRMatrix * A,HYPRE_Int * CF_marker,hypre_ParCSRMatrix * S,HYPRE_Int num_functions,HYPRE_Int * dof_func,HYPRE_Int skip_fine_or_same_sign)559 HYPRE_Int hypre_exchange_interp_data(
560     HYPRE_Int **CF_marker_offd,
561     HYPRE_Int **dof_func_offd,
562     hypre_CSRMatrix **A_ext,
563     HYPRE_Int *full_off_procNodes,
564     hypre_CSRMatrix **Sop,
565     hypre_ParCSRCommPkg **extend_comm_pkg,
566     hypre_ParCSRMatrix *A,
567     HYPRE_Int *CF_marker,
568     hypre_ParCSRMatrix *S,
569     HYPRE_Int num_functions,
570     HYPRE_Int *dof_func,
571     HYPRE_Int skip_fine_or_same_sign) // skip_fine_or_same_sign if we want to skip fine points in S and nnz with the same sign as diagonal in A
572 {
573 #ifdef HYPRE_PROFILE
574   hypre_profile_times[HYPRE_TIMER_ID_EXCHANGE_INTERP_DATA] -= hypre_MPI_Wtime();
575 #endif
576 
577   hypre_ParCSRCommPkg   *comm_pkg = hypre_ParCSRMatrixCommPkg(A);
578   hypre_CSRMatrix *A_diag = hypre_ParCSRMatrixDiag(A);
579   hypre_CSRMatrix *A_offd = hypre_ParCSRMatrixOffd(A);
580   HYPRE_Int              num_cols_A_offd = hypre_CSRMatrixNumCols(A_offd);
581   HYPRE_BigInt          *col_map_offd = hypre_ParCSRMatrixColMapOffd(A);
582   HYPRE_BigInt           col_1 = hypre_ParCSRMatrixFirstRowIndex(A);
583   HYPRE_Int              local_numrows = hypre_CSRMatrixNumRows(A_diag);
584   HYPRE_BigInt           col_n = col_1 + (HYPRE_BigInt)local_numrows;
585   HYPRE_BigInt          *found = NULL;
586 
587   /*----------------------------------------------------------------------
588    * Get the off processors rows for A and S, associated with columns in
589    * A_offd and S_offd.
590    *---------------------------------------------------------------------*/
591   *CF_marker_offd = hypre_TAlloc(HYPRE_Int, num_cols_A_offd, HYPRE_MEMORY_HOST);
592   hypre_exchange_marker(comm_pkg, CF_marker, *CF_marker_offd);
593 
594   hypre_ParCSRCommHandle *comm_handle_a_idx, *comm_handle_a_data;
595   *A_ext         = hypre_ParCSRMatrixExtractBExt_Overlap(A,A,1,&comm_handle_a_idx,&comm_handle_a_data,CF_marker,*CF_marker_offd,skip_fine_or_same_sign,skip_fine_or_same_sign);
596   HYPRE_Int *A_ext_i        = hypre_CSRMatrixI(*A_ext);
597   HYPRE_BigInt *A_ext_j        = hypre_CSRMatrixBigJ(*A_ext);
598   HYPRE_Int  A_ext_rows     = hypre_CSRMatrixNumRows(*A_ext);
599 
600   hypre_ParCSRCommHandle *comm_handle_s_idx;
601   *Sop           = hypre_ParCSRMatrixExtractBExt_Overlap(S,A,0,&comm_handle_s_idx,NULL,CF_marker,*CF_marker_offd,skip_fine_or_same_sign,0);
602   HYPRE_Int *Sop_i          = hypre_CSRMatrixI(*Sop);
603   HYPRE_BigInt *Sop_j       = hypre_CSRMatrixBigJ(*Sop);
604   HYPRE_Int  Soprows        = hypre_CSRMatrixNumRows(*Sop);
605 
606   HYPRE_Int *send_idx = (HYPRE_Int *)comm_handle_s_idx->send_data;
607   hypre_ParCSRCommHandleDestroy(comm_handle_s_idx);
608   hypre_TFree(send_idx, HYPRE_MEMORY_HOST);
609 
610   send_idx = (HYPRE_Int *)comm_handle_a_idx->send_data;
611   hypre_ParCSRCommHandleDestroy(comm_handle_a_idx);
612   hypre_TFree(send_idx, HYPRE_MEMORY_HOST);
613 
614   /* Find nodes that are neighbors of neighbors, not found in offd */
615 #ifdef HYPRE_PROFILE
616   hypre_profile_times[HYPRE_TIMER_ID_EXCHANGE_INTERP_DATA] += hypre_MPI_Wtime();
617 #endif
618   HYPRE_Int newoff = hypre_new_offd_nodes(&found, A_ext_rows, A_ext_i, A_ext_j,
619       Soprows, col_map_offd, col_1, col_n,
620       Sop_i, Sop_j, *CF_marker_offd);
621 #ifdef HYPRE_PROFILE
622   hypre_profile_times[HYPRE_TIMER_ID_EXCHANGE_INTERP_DATA] -= hypre_MPI_Wtime();
623 #endif
624   if(newoff >= 0)
625     *full_off_procNodes = newoff + num_cols_A_offd;
626   else
627   {
628     return hypre_error_flag;
629   }
630 
631   /* Possibly add new points and new processors to the comm_pkg, all
632    * processors need new_comm_pkg */
633 
634   /* AHB - create a new comm package just for extended info -
635      this will work better with the assumed partition*/
636   hypre_ParCSRFindExtendCommPkg(hypre_ParCSRMatrixComm(A),
637                                 hypre_ParCSRMatrixGlobalNumCols(A),
638                                 hypre_ParCSRMatrixFirstColDiag(A),
639                                 hypre_CSRMatrixNumCols(A_diag),
640                                 hypre_ParCSRMatrixColStarts(A),
641                                 hypre_ParCSRMatrixAssumedPartition(A),
642                                 newoff,
643                                 found,
644                                 extend_comm_pkg);
645 
646   *CF_marker_offd = hypre_TReAlloc(*CF_marker_offd, HYPRE_Int, *full_off_procNodes, HYPRE_MEMORY_HOST);
647   hypre_exchange_marker(*extend_comm_pkg, CF_marker, *CF_marker_offd + A_ext_rows);
648 
649   if(num_functions > 1)
650   {
651     if (*full_off_procNodes > 0)
652       *dof_func_offd = hypre_CTAlloc(HYPRE_Int, *full_off_procNodes, HYPRE_MEMORY_HOST);
653 
654     hypre_alt_insert_new_nodes(comm_pkg, *extend_comm_pkg, dof_func,
655         *full_off_procNodes, *dof_func_offd);
656   }
657 
658   hypre_TFree(found, HYPRE_MEMORY_HOST);
659 
660   HYPRE_Real *send_data = (HYPRE_Real *)comm_handle_a_data->send_data;
661   hypre_ParCSRCommHandleDestroy(comm_handle_a_data);
662   hypre_TFree(send_data, HYPRE_MEMORY_HOST);
663 
664 #ifdef HYPRE_PROFILE
665   hypre_profile_times[HYPRE_TIMER_ID_EXCHANGE_INTERP_DATA] += hypre_MPI_Wtime();
666 #endif
667 
668   return hypre_error_flag;
669 }
670 
hypre_build_interp_colmap(hypre_ParCSRMatrix * P,HYPRE_Int full_off_procNodes,HYPRE_Int * tmp_CF_marker_offd,HYPRE_BigInt * fine_to_coarse_offd)671 void hypre_build_interp_colmap(hypre_ParCSRMatrix *P, HYPRE_Int full_off_procNodes, HYPRE_Int *tmp_CF_marker_offd, HYPRE_BigInt *fine_to_coarse_offd)
672 {
673 #ifdef HYPRE_PROFILE
674    hypre_profile_times[HYPRE_TIMER_ID_RENUMBER_COLIDX] -= hypre_MPI_Wtime();
675 #endif
676    HYPRE_Int     n_fine = hypre_CSRMatrixNumRows(P->diag);
677 
678    HYPRE_Int     P_offd_size = P->offd->i[n_fine];
679    HYPRE_Int    *P_offd_j = P->offd->j;
680    HYPRE_BigInt *col_map_offd_P = NULL;
681    HYPRE_Int    *P_marker = NULL;
682    HYPRE_Int    *prefix_sum_workspace;
683    HYPRE_Int     num_cols_P_offd = 0;
684    HYPRE_Int     i, index;
685 
686    if (full_off_procNodes)
687    {
688       P_marker = hypre_TAlloc(HYPRE_Int, full_off_procNodes, HYPRE_MEMORY_HOST);
689    }
690    prefix_sum_workspace = hypre_TAlloc(HYPRE_Int, hypre_NumThreads() + 1, HYPRE_MEMORY_HOST);
691 
692 #ifdef HYPRE_USING_OPENMP
693 #pragma omp parallel for private(i) HYPRE_SMP_SCHEDULE
694 #endif
695    for (i = 0; i < full_off_procNodes; i++)
696    {
697       P_marker[i] = 0;
698    }
699 
700    /* These two loops set P_marker[i] to 1 if it appears in P_offd_j and if
701     * tmp_CF_marker_offd has i marked. num_cols_P_offd is then set to the
702     * total number of times P_marker is set */
703 #ifdef HYPRE_USING_OPENMP
704 #pragma omp parallel for private(i,index) HYPRE_SMP_SCHEDULE
705 #endif
706    for (i = 0; i < P_offd_size; i++)
707    {
708       index = P_offd_j[i];
709       if (tmp_CF_marker_offd[index] >= 0)
710       {
711          P_marker[index] = 1;
712       }
713    }
714 
715 #ifdef HYPRE_USING_OPENMP
716 #pragma omp parallel private(i)
717 #endif
718    {
719       HYPRE_Int i_begin, i_end;
720       hypre_GetSimpleThreadPartition(&i_begin, &i_end, full_off_procNodes);
721 
722       HYPRE_Int local_num_cols_P_offd = 0;
723       for (i = i_begin; i < i_end; i++)
724       {
725          if (P_marker[i] == 1) local_num_cols_P_offd++;
726       }
727 
728       hypre_prefix_sum(&local_num_cols_P_offd, &num_cols_P_offd, prefix_sum_workspace);
729 
730 #ifdef HYPRE_USING_OPENMP
731 #pragma omp master
732 #endif
733       {
734          if (num_cols_P_offd)
735          {
736             col_map_offd_P = hypre_TAlloc(HYPRE_BigInt, num_cols_P_offd, HYPRE_MEMORY_HOST);
737          }
738       }
739 #ifdef HYPRE_USING_OPENMP
740 #pragma omp barrier
741 #endif
742 
743       for (i = i_begin; i < i_end; i++)
744       {
745          if (P_marker[i] == 1)
746          {
747            col_map_offd_P[local_num_cols_P_offd++] = fine_to_coarse_offd[i];
748          }
749       }
750    }
751 
752    hypre_UnorderedBigIntMap col_map_offd_P_inverse;
753    hypre_big_sort_and_create_inverse_map(col_map_offd_P, num_cols_P_offd, &col_map_offd_P, &col_map_offd_P_inverse);
754 
755    // find old idx -> new idx map
756 #ifdef HYPRE_USING_OPENMP
757 #pragma omp parallel for
758 #endif
759    for (i = 0; i < full_off_procNodes; i++)
760    {
761       P_marker[i] = hypre_UnorderedBigIntMapGet(&col_map_offd_P_inverse, fine_to_coarse_offd[i]);
762    }
763 
764    if (num_cols_P_offd)
765    {
766       hypre_UnorderedBigIntMapDestroy(&col_map_offd_P_inverse);
767    }
768 
769 #ifdef HYPRE_USING_OPENMP
770 #pragma omp parallel for
771 #endif
772    for(i = 0; i < P_offd_size; i++)
773    {
774       P_offd_j[i] = P_marker[P_offd_j[i]];
775    }
776 
777    hypre_TFree(P_marker, HYPRE_MEMORY_HOST);
778    hypre_TFree(prefix_sum_workspace, HYPRE_MEMORY_HOST);
779 
780    if (num_cols_P_offd)
781    {
782       hypre_ParCSRMatrixColMapOffd(P) = col_map_offd_P;
783       hypre_CSRMatrixNumCols(P->offd) = num_cols_P_offd;
784    }
785 
786 #ifdef HYPRE_PROFILE
787    hypre_profile_times[HYPRE_TIMER_ID_RENUMBER_COLIDX] += hypre_MPI_Wtime();
788 #endif
789 }
790