← Back to branch summary

~ubuntu-branches/ubuntu/precise/code-saturne/precise

« back to all changes in this revision

Viewing changes to src/base/cs_perio.c

  • Committer: Package Import Robot
  • Author(s): Sylvestre Ledru
  • Date: 2011-11-01 17:43:32 UTC
  • mto: (6.1.7 sid)
  • mto: This revision was merged to the branch mainline in revision 11.
  • Revision ID: package-import@ubuntu.com-20111101174332-tl4vk45no0x3emc3
Tags: upstream-2.1.0
Import upstream version 2.1.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/base/cs_perio.c
1
1
/*============================================================================
2
 
 *
3
 
 *     This file is part of the Code_Saturne Kernel, element of the
4
 
 *     Code_Saturne CFD tool.
5
 
 *
6
 
 *     Copyright (C) 1998-2009 EDF S.A., France
7
 
 *
8
 
 *     contact: saturne-support@edf.fr
9
 
 *
10
 
 *     The Code_Saturne Kernel is free software; you can redistribute it
11
 
 *     and/or modify it under the terms of the GNU General Public License
12
 
 *     as published by the Free Software Foundation; either version 2 of
13
 
 *     the License, or (at your option) any later version.
14
 
 *
15
 
 *     The Code_Saturne Kernel is distributed in the hope that it will be
16
 
 *     useful, but WITHOUT ANY WARRANTY; without even the implied warranty
17
 
 *     of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18
 
 *     GNU General Public License for more details.
19
 
 *
20
 
 *     You should have received a copy of the GNU General Public License
21
 
 *     along with the Code_Saturne Kernel; if not, write to the
22
 
 *     Free Software Foundation, Inc.,
23
 
 *     51 Franklin St, Fifth Floor,
24
 
 *     Boston, MA  02110-1301  USA
25
 
 *
26
 
 *============================================================================*/
27
 
 
28
 
/*============================================================================
29
 
 * Functions handling with periodicity.
30
 
 *============================================================================*/
 
2
 * Functions handling periodicity.
 
3
 *============================================================================*/
 
4
 
 
5
/*
 
6
  This file is part of Code_Saturne, a general-purpose CFD tool.
 
7
 
 
8
  Copyright (C) 1998-2011 EDF S.A.
 
9
 
 
10
  This program is free software; you can redistribute it and/or modify it under
 
11
  the terms of the GNU General Public License as published by the Free Software
 
12
  Foundation; either version 2 of the License, or (at your option) any later
 
13
  version.
 
14
 
 
15
  This program is distributed in the hope that it will be useful, but WITHOUT
 
16
  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 
17
  FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
 
18
  details.
 
19
 
 
20
  You should have received a copy of the GNU General Public License along with
 
21
  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
 
22
  Street, Fifth Floor, Boston, MA 02110-1301, USA.
 
23
*/
 
24
 
 
25
/*----------------------------------------------------------------------------*/
31
26
 
32
27
#if defined(HAVE_CONFIG_H)
33
28
#include "cs_config.h"
56
51
 * FVM library headers
57
52
 *----------------------------------------------------------------------------*/
58
53
 
 
54
#include <fvm_interface.h>
59
55
#include <fvm_periodicity.h>
60
56
#include <fvm_parall.h>
61
57
 
81
77
 * Local Macro definitions
82
78
 *============================================================================*/
83
79
 
 
80
/*=============================================================================
 
81
 * Local Structure definitions
 
82
 *============================================================================*/
 
83
 
 
84
/* Structure used for building mesh structure */
 
85
/* ------------------------------------------ */
 
86
 
 
87
typedef struct {
 
88
 
 
89
  /* Periodic features */
 
90
 
 
91
  fvm_lnum_t  *per_face_idx;    /* Index on periodicity for per_face_lst */
 
92
 
 
93
  fvm_lnum_t  *per_face_lst;    /* Periodic faces list. For each couple,
 
94
                                   we have the local face number on local rank
 
95
                                   and the local face number on distant rank */
 
96
 
 
97
  cs_int_t   *per_rank_lst;     /* Remote ranks list. For each couple,
 
98
                                   we have the distant rank number. Exists
 
99
                                   only in case of parallelism. */
 
100
 
 
101
} _perio_mesh_builder_t ;
 
102
 
84
103
/*============================================================================
85
104
 * Static global variables
86
105
 *============================================================================*/
251
270
}
252
271
 
253
272
/*----------------------------------------------------------------------------
 
273
 * Compute a matrix/vector product to apply a transformation to a given
 
274
 * vector.
 
275
 *
 
276
 * parameters:
 
277
 *   matrix[3][4] --> matrix of the transformation in homogeneous coord.
 
278
 *                    last line = [0; 0; 0; 1] (Not used here)
 
279
 *   xyz          <-> array of coordinates
 
280
 *----------------------------------------------------------------------------*/
 
281
 
 
282
static void
 
283
_apply_vector_rotation_i(cs_real_t    matrix[3][4],
 
284
                         cs_real_t   *xyz)
 
285
{
 
286
  cs_int_t  i;
 
287
 
 
288
  cs_real_t  t[3];
 
289
  for (i = 0; i < 3; i++)
 
290
    t[i] = xyz[i];
 
291
 
 
292
  /* Initialize output */
 
293
 
 
294
  for (i = 0; i < 3; i++)
 
295
    xyz[i] = matrix[i][0]*t[0] + matrix[i][1]*t[1] + matrix[i][2]*t[2];
 
296
 
 
297
}
 
298
 
 
299
/*----------------------------------------------------------------------------
254
300
 * Compute a matrix * tensor * Tmatrix product to apply a rotation to a
255
301
 * given tensor
256
302
 *
356
402
 *   cell_id   --> cell id
357
403
 *   rom       --> density array
358
404
 *   call_id   --> first or second call
359
 
 *   phase_id  --> phase id
360
405
 *   dudxyz    <-> gradient on the components of the velocity.
361
406
 *   wdudxy    <-> associated working array.
362
407
 *----------------------------------------------------------------------------*/
366
411
               cs_int_t          cell_id,
367
412
               const cs_real_t  *rom,
368
413
               cs_int_t          call_id,
369
 
               cs_int_t          phase_id,
370
414
               cs_real_t        *dudxyz,
371
415
               cs_real_t        *wdudxy)
372
416
{
381
425
 
382
426
    for (i = 0; i < 3; i++) {
383
427
      for (j = 0; j < 3; j++) {
384
 
        id = h_cell_id + n_ghost_cells*i + stride*j + 3*stride*phase_id;
 
428
        id = h_cell_id + n_ghost_cells*i + stride*j;
385
429
        wdudxy[id] = dudxyz[id];
386
430
        dudxyz[id] *= rom[cell_id];
387
431
      }
392
436
 
393
437
    for (i = 0; i < 3; i++) {
394
438
      for (j = 0; j < 3; j++) {
395
 
        id = h_cell_id + n_ghost_cells*i + stride*j + 3*stride*phase_id;
 
439
        id = h_cell_id + n_ghost_cells*i + stride*j;
396
440
        dudxyz[id] = wdudxy[id];
397
441
      }
398
442
    }
410
454
 *   cell_id   --> cell id
411
455
 *   rom       --> density array
412
456
 *   call_id   --> first or second call
413
 
 *   phase_id  --> phase id
414
457
 *   drdxyz    <-> Gradient on components of Rij (Reynolds stress tensor)
415
458
 *   wdrdxy    <-> associated working array.
416
459
 *----------------------------------------------------------------------------*/
420
463
               cs_int_t          cell_id,
421
464
               const cs_real_t  *rom,
422
465
               cs_int_t          call_id,
423
 
               cs_int_t          phase_id,
424
466
               cs_real_t        *drdxyz,
425
467
               cs_real_t        *wdrdxy)
426
468
{
435
477
 
436
478
    for (i = 0; i < 2*3; i++) {
437
479
      for (j = 0; j < 3; j++) {
438
 
        id = h_cell_id + n_ghost_cells*i + stride*j +3*stride*phase_id;
 
480
        id = h_cell_id + n_ghost_cells*i + stride*j;
439
481
        wdrdxy[id] = drdxyz[id];
440
482
        drdxyz[id] *= rom[cell_id];
441
483
      }
446
488
 
447
489
    for (i = 0; i < 2*3; i++) {
448
490
      for (j = 0; j < 3; j++) {
449
 
        id = h_cell_id + n_ghost_cells*i + stride*j +3*stride*phase_id;
 
491
        id = h_cell_id + n_ghost_cells*i + stride*j;
450
492
        drdxyz[id] = wdrdxy[id];
451
493
      }
452
494
    }
460
502
 * parameters:
461
503
 *   strid_c    --> stride on the component
462
504
 *   strid_v    --> stride on the variable
463
 
 *   strid_p    --> stride on the phase
464
505
 *   dxyz       <-> gradient on the variable (dudxy or drdxy)
465
506
 *   w1, w2, w3 <-> working buffers
466
507
 *----------------------------------------------------------------------------*/
468
509
static void
469
510
_peinur1(cs_int_t      strid_c,
470
511
         cs_int_t      strid_v,
471
 
         cs_int_t      strid_p,
472
512
         cs_real_t    *dxyz,
473
513
         cs_real_t    *w1,
474
514
         cs_real_t    *w2,
475
515
         cs_real_t    *w3)
476
516
{
477
517
  cs_int_t  i, t_id, rank_id, shift;
478
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
518
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
479
519
 
480
520
  cs_mesh_t  *mesh = cs_glob_mesh;
481
521
  cs_halo_t  *halo = mesh->halo;
515
555
      }
516
556
 
517
557
      for (i = start_std; i < end_std; i++) {
518
 
        dxyz[i + strid_c + strid_v*0 + strid_p] = w1[n_cells + i];
519
 
        dxyz[i + strid_c + strid_v*1 + strid_p] = w2[n_cells + i];
520
 
        dxyz[i + strid_c + strid_v*2 + strid_p] = w3[n_cells + i];
 
558
        dxyz[i + strid_c + strid_v*0] = w1[n_cells + i];
 
559
        dxyz[i + strid_c + strid_v*1] = w2[n_cells + i];
 
560
        dxyz[i + strid_c + strid_v*2] = w3[n_cells + i];
521
561
      }
522
562
 
523
563
      if (mesh->halo_type == CS_HALO_EXTENDED) {
524
564
 
525
565
        for (i = start_ext; i < end_ext; i++) {
526
 
          dxyz[i + strid_c + strid_v*0 + strid_p] = w1[n_cells + i];
527
 
          dxyz[i + strid_c + strid_v*1 + strid_p] = w2[n_cells + i];
528
 
          dxyz[i + strid_c + strid_v*2 + strid_p] = w3[n_cells + i];
 
566
          dxyz[i + strid_c + strid_v*0] = w1[n_cells + i];
 
567
          dxyz[i + strid_c + strid_v*1] = w2[n_cells + i];
 
568
          dxyz[i + strid_c + strid_v*2] = w3[n_cells + i];
529
569
        }
530
570
 
531
571
      } /* End if extended halo */
636
676
 
637
677
  const cs_halo_t *halo = cs_glob_mesh->halo;
638
678
 
 
679
  if (halo == NULL)
 
680
    return;
 
681
 
639
682
  /* 1. Checking    */
640
683
  /*----------------*/
641
684
 
687
730
       rotations, something has been done before; see PERINR for example). */
688
731
 
689
732
    else if (*itenso == 2)
690
 
      cs_perio_sync_var_vect(halo,
691
 
                             CS_HALO_STANDARD,
692
 
                             CS_PERIO_ROTA_IGNORE,
693
 
                             var11, var22, var33);
 
733
      cs_perio_sync_var_vect_ni(halo,
 
734
                                CS_HALO_STANDARD,
 
735
                                CS_PERIO_ROTA_IGNORE,
 
736
                                var11, var22, var33);
694
737
 
695
738
  } /* End of idimte == 0 case */
696
739
 
698
741
     it is (at least) a vector. Translation and rotation are exchanged. */
699
742
 
700
743
  else if (*idimte == 1)
701
 
    cs_perio_sync_var_vect(halo,
702
 
                           CS_HALO_STANDARD,
703
 
                           CS_PERIO_ROTA_COPY,
704
 
                           var11, var22, var33);
 
744
    cs_perio_sync_var_vect_ni(halo,
 
745
                              CS_HALO_STANDARD,
 
746
                              CS_PERIO_ROTA_COPY,
 
747
                              var11, var22, var33);
705
748
 
706
749
  /* --> If we want to handle the variable as a tensor, we suppose that
707
750
     it is a tensor. Translation and rotation are exchanged. */
708
751
 
709
752
  else if (*idimte == 2)
710
 
    cs_perio_sync_var_tens(halo,
711
 
                           CS_HALO_STANDARD,
712
 
                           var11, var12, var13,
713
 
                           var21, var22, var23,
714
 
                           var31, var32, var33);
 
753
    cs_perio_sync_var_tens_ni(halo,
 
754
                              CS_HALO_STANDARD,
 
755
                              var11, var12, var13,
 
756
                              var21, var22, var23,
 
757
                              var31, var32, var33);
715
758
 
716
759
  /* --> If we want to handle the variable as a tensor, but that
717
760
     it is a tensor's diagonal, we suppose that it is a tensor.
718
761
     Translation and rotation are exchanged. */
719
762
 
720
763
  else if (*idimte == 21)
721
 
    cs_perio_sync_var_diag(halo,
722
 
                           CS_HALO_STANDARD,
723
 
                           var11, var22, var33);
 
764
    cs_perio_sync_var_diag_ni(halo,
 
765
                              CS_HALO_STANDARD,
 
766
                              var11, var22, var33);
724
767
}
725
768
 
726
769
/*----------------------------------------------------------------------------
829
872
       rotations, something has been done before; see PERINR for example). */
830
873
 
831
874
    else if (*itenso == 2)
832
 
      cs_perio_sync_var_vect(halo,
833
 
                             CS_HALO_EXTENDED,
834
 
                             CS_PERIO_ROTA_IGNORE,
835
 
                             var11, var22, var33);
 
875
      cs_perio_sync_var_vect_ni(halo,
 
876
                                CS_HALO_EXTENDED,
 
877
                                CS_PERIO_ROTA_IGNORE,
 
878
                                var11, var22, var33);
836
879
 
837
880
  } /* End of idimte == 0 case */
838
881
 
840
883
     it is (at least) a vector. Translation and rotation are exchanged. */
841
884
 
842
885
  else if (*idimte == 1)
843
 
    cs_perio_sync_var_vect(halo,
844
 
                           CS_HALO_EXTENDED,
845
 
                           CS_PERIO_ROTA_COPY,
846
 
                           var11, var22, var33);
 
886
    cs_perio_sync_var_vect_ni(halo,
 
887
                              CS_HALO_EXTENDED,
 
888
                              CS_PERIO_ROTA_COPY,
 
889
                              var11, var22, var33);
847
890
 
848
891
  /* --> If we want to handle the variable as a tensor, we suppose that
849
892
     it is a tensor. Translation and rotation are exchanged. */
850
893
 
851
894
  else if (*idimte == 2)
852
 
    cs_perio_sync_var_tens(halo,
853
 
                           CS_HALO_EXTENDED,
854
 
                           var11, var12, var13,
855
 
                           var21, var22, var23,
856
 
                           var31, var32, var33);
 
895
    cs_perio_sync_var_tens_ni(halo,
 
896
                              CS_HALO_EXTENDED,
 
897
                              var11, var12, var13,
 
898
                              var21, var22, var23,
 
899
                              var31, var32, var33);
857
900
 
858
901
  /* --> If we want to handle the variable as a tensor, but that
859
902
     it is a tensor's diagonal, we suppose that it is a tensor.
860
903
     Translation and rotation are exchanged. */
861
904
 
862
905
  else if (*idimte == 21)
863
 
    cs_perio_sync_var_diag(halo,
864
 
                           CS_HALO_EXTENDED,
865
 
                           var11, var22, var33);
 
906
    cs_perio_sync_var_diag_ni(halo,
 
907
                              CS_HALO_EXTENDED,
 
908
                              var11, var22, var33);
866
909
 
867
910
}
868
911
 
885
928
 * INTEGER          IMASPE      :  -> : suivant l'appel de INIMAS
886
929
 *                                          = 1 si appel de RESOLP ou NAVSTO
887
930
 *                                          = 2 si appel de DIVRIJ
888
 
 * INTEGER          IPHAS       :  -> : numero de phase courante
889
931
 * INTEGER          IMASPE      :  -> : indicateur d'appel dans INIMAS
890
932
 *                                          = 1 si appel au debut
891
933
 *                                          = 2 si appel a la fin
897
939
 * DOUBLE PRECISION WDUDXY      :  -  : tableau de travail pour DUDXYZ
898
940
 * DOUBLE PRECISION WDRDXY      :  -  : tableau de travail pour DRDXYZ
899
941
 *
900
 
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3*NPHAS
901
 
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3*NPHAS
 
942
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3
 
943
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3
902
944
 *----------------------------------------------------------------------------*/
903
945
 
904
946
void
905
947
CS_PROCF (permas, PERMAS)(const cs_int_t    *imaspe,
906
 
                          const cs_int_t    *iphas,
907
948
                          const cs_int_t    *iappel,
908
949
                          cs_real_t          rom[],
909
950
                          cs_real_t         *dudxyz,
912
953
                          cs_real_t         *wdrdxy)
913
954
{
914
955
  cs_int_t  i, cell_id, rank_id, shift, t_id;
915
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
956
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
916
957
 
917
958
  cs_mesh_t  *mesh = cs_glob_mesh;
918
959
  cs_halo_t  *halo = mesh->halo;
919
960
  cs_halo_type_t  halo_type = mesh->halo_type;
920
961
 
921
 
  const cs_int_t  phase_id = *iphas - 1;
922
 
 
923
962
  if (halo_type == CS_HALO_N_TYPES)
924
963
    return;
925
964
 
949
988
        cell_id = mesh->n_cells + i;
950
989
 
951
990
        if (*imaspe == 1)
952
 
          _update_dudxyz(i, cell_id, rom, *iappel, phase_id, dudxyz, wdudxy);
 
991
          _update_dudxyz(i, cell_id, rom, *iappel, dudxyz, wdudxy);
953
992
 
954
993
        if (*imaspe == 2)
955
 
          _update_drdxyz(i, cell_id, rom, *iappel, phase_id, drdxyz, wdrdxy);
 
994
          _update_drdxyz(i, cell_id, rom, *iappel, drdxyz, wdrdxy);
956
995
 
957
996
      } /* End of loop on halo elements */
958
997
 
963
1002
          cell_id = mesh->n_cells + i;
964
1003
 
965
1004
          if (*imaspe == 1)
966
 
            _update_dudxyz(i, cell_id, rom, *iappel, phase_id, dudxyz, wdudxy);
 
1005
            _update_dudxyz(i, cell_id, rom, *iappel, dudxyz, wdudxy);
967
1006
 
968
1007
          if (*imaspe == 2)
969
 
            _update_drdxyz(i, cell_id, rom, *iappel, phase_id, drdxyz, wdrdxy);
 
1008
            _update_drdxyz(i, cell_id, rom, *iappel, drdxyz, wdrdxy);
970
1009
 
971
1010
        } /* End of loop on halo elements */
972
1011
 
1003
1042
 * SUBROUTINE PERING
1004
1043
 * *****************
1005
1044
 *
1006
 
 * INTEGER          NPHAS        :  -> : numero de phase courante
1007
1045
 * INTEGER          IVAR         :  -> : numero de la variable
1008
1046
 * INTEGER          IDIMTE       : <-  : dimension de la variable (maximum 3)
1009
1047
 *                                        0 : scalaire (VAR11), ou assimile
1046
1084
 * DOUBLE PRECISION DRDXYZ       :  -> : gradient de R aux cellules halo pour
1047
1085
 *                                       l'approche explicite en periodicite
1048
1086
 *
1049
 
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3*NPHAS
1050
 
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3*NPHAS
 
1087
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3
 
1088
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3
1051
1089
 *----------------------------------------------------------------------------*/
1052
1090
 
1053
1091
void
1054
 
CS_PROCF (pering, PERING)(const cs_int_t    *nphas,
1055
 
                          const cs_int_t    *ivar,
 
1092
CS_PROCF (pering, PERING)(const cs_int_t    *ivar,
1056
1093
                          cs_int_t          *idimte,
1057
1094
                          cs_int_t          *itenso,
1058
1095
                          const cs_int_t    *iperot,
1059
1096
                          const cs_int_t    *iguper,
1060
1097
                          const cs_int_t    *igrper,
1061
 
                          const cs_int_t     iu[CS_NPHSMX],
1062
 
                          const cs_int_t     iv[CS_NPHSMX],
1063
 
                          const cs_int_t     iw[CS_NPHSMX],
1064
 
                          const cs_int_t     itytur[CS_NPHSMX],
1065
 
                          const cs_int_t     ir11[CS_NPHSMX],
1066
 
                          const cs_int_t     ir22[CS_NPHSMX],
1067
 
                          const cs_int_t     ir33[CS_NPHSMX],
1068
 
                          const cs_int_t     ir12[CS_NPHSMX],
1069
 
                          const cs_int_t     ir13[CS_NPHSMX],
1070
 
                          const cs_int_t     ir23[CS_NPHSMX],
 
1098
                          const cs_int_t    *iu,
 
1099
                          const cs_int_t    *iv,
 
1100
                          const cs_int_t    *iw,
 
1101
                          const cs_int_t    *itytur,
 
1102
                          const cs_int_t    *ir11,
 
1103
                          const cs_int_t    *ir22,
 
1104
                          const cs_int_t    *ir33,
 
1105
                          const cs_int_t    *ir12,
 
1106
                          const cs_int_t    *ir13,
 
1107
                          const cs_int_t    *ir23,
1071
1108
                          cs_real_t          dpdx[],
1072
1109
                          cs_real_t          dpdy[],
1073
1110
                          cs_real_t          dpdz[],
1074
1111
                          const cs_real_t   *dudxyz,
1075
1112
                          const cs_real_t   *drdxyz)
1076
1113
{
1077
 
  cs_int_t  i, rank_id, phase_id, t_id, shift;
1078
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext, tag;
1079
 
  cs_int_t  d_ph, d_var;
 
1114
  cs_int_t  i, rank_id, t_id, shift, tag;
 
1115
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
 
1116
  cs_int_t  d_ph = 0, d_var = 0;
1080
1117
 
1081
1118
  cs_mesh_t  *mesh = cs_glob_mesh;
1082
1119
  cs_halo_t  *halo = mesh->halo;
1085
1122
  const cs_int_t  n_transforms = mesh->n_transforms;
1086
1123
  const cs_int_t  n_ghost_cells = mesh->n_ghost_cells;
1087
1124
  const cs_int_t  stride1 = n_ghost_cells * 3;
1088
 
  const cs_int_t  stride2 = stride1 * 3;
 
1125
 
 
1126
  if (halo == NULL)
 
1127
    return;
1089
1128
 
1090
1129
  /* We treat the gradient like a vector by default ...
1091
1130
     (i.e. we assume that this is the gradient of a scalar. */
1109
1148
 
1110
1149
    tag = 0;
1111
1150
 
1112
 
    for (phase_id = 0; phase_id < *nphas; phase_id++) {
1113
 
 
1114
 
      if (   *ivar == iu[phase_id]
1115
 
          || *ivar == iv[phase_id]
1116
 
          || *ivar == iw[phase_id]) {
1117
 
 
1118
 
        d_ph = stride2*phase_id;
1119
 
        tag = 1;
1120
 
 
1121
 
        if (*ivar == iu[phase_id]) d_var = 0;
1122
 
        if (*ivar == iv[phase_id]) d_var = n_ghost_cells;
1123
 
        if (*ivar == iw[phase_id]) d_var = 2*n_ghost_cells;
1124
 
 
1125
 
        if (*iguper == 1) { /* dudxyz not computed */
1126
 
 
1127
 
          for (t_id = 0; t_id < n_transforms; t_id++) {
1128
 
 
1129
 
            shift = 4 * halo->n_c_domains * t_id;
1130
 
 
1131
 
            for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
1132
 
 
1133
 
              start_std = halo->perio_lst[shift + 4*rank_id];
1134
 
              length = halo->perio_lst[shift + 4*rank_id + 1];
1135
 
              end_std = start_std + length;
1136
 
 
1137
 
              if (mesh->halo_type == CS_HALO_EXTENDED) {
1138
 
 
1139
 
                start_ext = halo->perio_lst[shift + 4*rank_id + 2];
1140
 
                length = halo->perio_lst[shift + 4*rank_id + 3];
1141
 
                end_ext = start_ext + length;
1142
 
 
1143
 
              }
1144
 
 
1145
 
              for (i = start_std; i < end_std; i++) {
1146
 
                dpdx[n_cells + i] = dudxyz[i + d_var + stride1*0 + d_ph];
1147
 
                dpdy[n_cells + i] = dudxyz[i + d_var + stride1*1 + d_ph];
1148
 
                dpdz[n_cells + i] = dudxyz[i + d_var + stride1*2 + d_ph];
1149
 
              }
1150
 
 
1151
 
              if (mesh->halo_type == CS_HALO_EXTENDED) {
1152
 
 
1153
 
                for (i = start_ext; i < end_ext; i++) {
1154
 
                  dpdx[n_cells + i] = dudxyz[i + d_var + stride1*0 + d_ph];
1155
 
                  dpdy[n_cells + i] = dudxyz[i + d_var + stride1*1 + d_ph];
1156
 
                  dpdz[n_cells + i] = dudxyz[i + d_var + stride1*2 + d_ph];
1157
 
                }
1158
 
 
1159
 
              } /* End if extended halo */
1160
 
 
1161
 
            } /* End of loop on ranks */
1162
 
 
1163
 
          } /* End of loop on transformations */
1164
 
 
1165
 
        } /* End if *iguper == 1 */
1166
 
 
1167
 
      } /* If *ivar == iu or iv or iw */
1168
 
 
1169
 
      else if ((itytur[phase_id] == 3) &&
1170
 
               (*ivar == ir11[phase_id] || *ivar == ir22[phase_id] ||
1171
 
                *ivar == ir33[phase_id] || *ivar == ir12[phase_id] ||
1172
 
                *ivar == ir13[phase_id] || *ivar == ir23[phase_id])) {
1173
 
 
1174
 
        d_ph = 2*stride2*phase_id;
1175
 
        tag = 1;
1176
 
 
1177
 
        if (*ivar == ir11[phase_id]) d_var = 0;
1178
 
        if (*ivar == ir22[phase_id]) d_var = n_ghost_cells;
1179
 
        if (*ivar == ir33[phase_id]) d_var = 2*n_ghost_cells;
1180
 
        if (*ivar == ir12[phase_id]) d_var = 3*n_ghost_cells;
1181
 
        if (*ivar == ir13[phase_id]) d_var = 4*n_ghost_cells;
1182
 
        if (*ivar == ir23[phase_id]) d_var = 5*n_ghost_cells;
1183
 
 
1184
 
        if (*igrper == 1) {
1185
 
 
1186
 
          for (t_id = 0; t_id < n_transforms; t_id++) {
1187
 
 
1188
 
            shift = 4 * halo->n_c_domains * t_id;
1189
 
 
1190
 
            for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
1191
 
 
1192
 
              start_std = halo->perio_lst[shift + 4*rank_id];
1193
 
              length = halo->perio_lst[shift + 4*rank_id + 1];
1194
 
              end_std = start_std + length;
1195
 
 
1196
 
              if (mesh->halo_type == CS_HALO_EXTENDED) {
1197
 
 
1198
 
                start_ext = halo->perio_lst[shift + 4*rank_id + 2];
1199
 
                length = halo->perio_lst[shift + 4*rank_id + 3];
1200
 
                end_ext = start_ext + length;
1201
 
 
1202
 
              }
1203
 
 
1204
 
              for (i = start_std; i < end_std; i++) {
 
1151
    if (*ivar == *iu || *ivar == *iv || *ivar == *iw) {
 
1152
 
 
1153
      tag = 1;
 
1154
 
 
1155
      if (*ivar == *iu) d_var = 0;
 
1156
      if (*ivar == *iv) d_var = n_ghost_cells;
 
1157
      if (*ivar == *iw) d_var = 2*n_ghost_cells;
 
1158
 
 
1159
      if (*iguper == 1) { /* dudxyz not computed */
 
1160
 
 
1161
        for (t_id = 0; t_id < n_transforms; t_id++) {
 
1162
 
 
1163
          shift = 4 * halo->n_c_domains * t_id;
 
1164
 
 
1165
          for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
 
1166
 
 
1167
            start_std = halo->perio_lst[shift + 4*rank_id];
 
1168
            length = halo->perio_lst[shift + 4*rank_id + 1];
 
1169
            end_std = start_std + length;
 
1170
 
 
1171
            if (mesh->halo_type == CS_HALO_EXTENDED) {
 
1172
 
 
1173
              start_ext = halo->perio_lst[shift + 4*rank_id + 2];
 
1174
              length = halo->perio_lst[shift + 4*rank_id + 3];
 
1175
              end_ext = start_ext + length;
 
1176
 
 
1177
            }
 
1178
 
 
1179
            for (i = start_std; i < end_std; i++) {
 
1180
              dpdx[n_cells + i] = dudxyz[i + d_var + stride1*0];
 
1181
              dpdy[n_cells + i] = dudxyz[i + d_var + stride1*1];
 
1182
              dpdz[n_cells + i] = dudxyz[i + d_var + stride1*2];
 
1183
            }
 
1184
 
 
1185
            if (mesh->halo_type == CS_HALO_EXTENDED) {
 
1186
 
 
1187
              for (i = start_ext; i < end_ext; i++) {
 
1188
                dpdx[n_cells + i] = dudxyz[i + d_var + stride1*0];
 
1189
                dpdy[n_cells + i] = dudxyz[i + d_var + stride1*1];
 
1190
                dpdz[n_cells + i] = dudxyz[i + d_var + stride1*2];
 
1191
              }
 
1192
 
 
1193
            } /* End if extended halo */
 
1194
 
 
1195
          } /* End of loop on ranks */
 
1196
 
 
1197
        } /* End of loop on transformations */
 
1198
 
 
1199
      } /* End if *iguper == 1 */
 
1200
 
 
1201
    } /* If *ivar == iu or iv or iw */
 
1202
 
 
1203
    else if ((*itytur == 3) &&
 
1204
             (*ivar == *ir11 || *ivar == *ir22 ||
 
1205
              *ivar == *ir33 || *ivar == *ir12 ||
 
1206
              *ivar == *ir13 || *ivar == *ir23)) {
 
1207
 
 
1208
      tag = 1;
 
1209
 
 
1210
      if (*ivar == *ir11) d_var = 0;
 
1211
      if (*ivar == *ir22) d_var = n_ghost_cells;
 
1212
      if (*ivar == *ir33) d_var = 2*n_ghost_cells;
 
1213
      if (*ivar == *ir12) d_var = 3*n_ghost_cells;
 
1214
      if (*ivar == *ir13) d_var = 4*n_ghost_cells;
 
1215
      if (*ivar == *ir23) d_var = 5*n_ghost_cells;
 
1216
 
 
1217
      if (*igrper == 1) {
 
1218
 
 
1219
        for (t_id = 0; t_id < n_transforms; t_id++) {
 
1220
 
 
1221
          shift = 4 * halo->n_c_domains * t_id;
 
1222
 
 
1223
          for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
 
1224
 
 
1225
            start_std = halo->perio_lst[shift + 4*rank_id];
 
1226
            length = halo->perio_lst[shift + 4*rank_id + 1];
 
1227
            end_std = start_std + length;
 
1228
 
 
1229
            if (mesh->halo_type == CS_HALO_EXTENDED) {
 
1230
 
 
1231
              start_ext = halo->perio_lst[shift + 4*rank_id + 2];
 
1232
              length = halo->perio_lst[shift + 4*rank_id + 3];
 
1233
              end_ext = start_ext + length;
 
1234
 
 
1235
            }
 
1236
 
 
1237
            for (i = start_std; i < end_std; i++) {
 
1238
              dpdx[n_cells + i] = drdxyz[i + d_var + 2*stride1*0 + d_ph];
 
1239
              dpdy[n_cells + i] = drdxyz[i + d_var + 2*stride1*1 + d_ph];
 
1240
              dpdz[n_cells + i] = drdxyz[i + d_var + 2*stride1*2 + d_ph];
 
1241
            }
 
1242
 
 
1243
            if (mesh->halo_type == CS_HALO_EXTENDED) {
 
1244
 
 
1245
              for (i = start_ext; i < end_ext; i++) {
1205
1246
                dpdx[n_cells + i] = drdxyz[i + d_var + 2*stride1*0 + d_ph];
1206
1247
                dpdy[n_cells + i] = drdxyz[i + d_var + 2*stride1*1 + d_ph];
1207
1248
                dpdz[n_cells + i] = drdxyz[i + d_var + 2*stride1*2 + d_ph];
1208
1249
              }
1209
1250
 
1210
 
              if (mesh->halo_type == CS_HALO_EXTENDED) {
1211
 
 
1212
 
                for (i = start_ext; i < end_ext; i++) {
1213
 
                  dpdx[n_cells + i] = drdxyz[i + d_var + 2*stride1*0 + d_ph];
1214
 
                  dpdy[n_cells + i] = drdxyz[i + d_var + 2*stride1*1 + d_ph];
1215
 
                  dpdz[n_cells + i] = drdxyz[i + d_var + 2*stride1*2 + d_ph];
1216
 
                }
1217
 
 
1218
 
              } /* End if extended halo */
1219
 
 
1220
 
            } /* End of loop on ranks */
1221
 
 
1222
 
          } /* End of loop on transformations */
1223
 
 
1224
 
        } /* End if *igrper == 1 */
1225
 
 
1226
 
      } /* If itytur[phase_id] == 3 and *ivar == irij */
1227
 
 
1228
 
    } /* End of loop on phases */
 
1251
            } /* End if extended halo */
 
1252
 
 
1253
          } /* End of loop on ranks */
 
1254
 
 
1255
        } /* End of loop on transformations */
 
1256
 
 
1257
      } /* End if *igrper == 1 */
 
1258
 
 
1259
    } /* If *itytur == 3 and *ivar == irij */
1229
1260
 
1230
1261
    if (tag == 1) {
1231
1262
      *idimte = 0;
1244
1275
 * *****************
1245
1276
 *
1246
1277
 * INTEGER          ISOU          :  -> : component of the velocity vector
1247
 
 * INTEGER          IPHAS         :  -> : current phase number
1248
1278
 * DOUBLE PRECISION DUDXYZ        : <-> : gradient of the velocity vector
1249
1279
 *                                        for ghost cells and for an explicit
1250
1280
 *                                        treatment of the periodicity.
1251
1281
 * DOUBLE PRECISION W1..3(NCELET) :  -  : working buffers
1252
1282
 *
1253
 
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3*NPHAS
 
1283
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3
1254
1284
 *----------------------------------------------------------------------------*/
1255
1285
 
1256
1286
void
1257
1287
CS_PROCF (peinu1, PEINU1)(const cs_int_t    *isou,
1258
 
                          const cs_int_t    *iphas,
1259
1288
                          cs_real_t         *dudxyz,
1260
1289
                          cs_real_t          w1[],
1261
1290
                          cs_real_t          w2[],
1265
1294
 
1266
1295
  const cs_int_t  n_ghost_cells = mesh->n_ghost_cells;
1267
1296
  const cs_int_t  comp_id = *isou - 1;
1268
 
  const cs_int_t  phas_id = *iphas - 1;
1269
1297
  const cs_int_t  strid_v = n_ghost_cells * 3;
1270
 
  const cs_int_t  strid_p = strid_v * 3 * phas_id;
1271
1298
  const cs_int_t  strid_c = n_ghost_cells * comp_id;
1272
1299
 
1273
 
  _peinur1(strid_c, strid_v, strid_p, dudxyz, w1, w2, w3);
 
1300
  _peinur1(strid_c, strid_v, dudxyz, w1, w2, w3);
1274
1301
}
1275
1302
 
1276
1303
/*----------------------------------------------------------------------------
1281
1308
 * SUBROUTINE PEINU2 (VAR)
1282
1309
 * *****************
1283
1310
 *
1284
 
 * INTEGER          IPHAS         :  -> : current phase number
1285
1311
 * DOUBLE PRECISION DUDXYZ        : <-> : gradient of the velocity vector
1286
1312
 *                                        for ghost cells and for an explicit
1287
1313
 *                                        treatment of the periodicity.
1288
1314
 *
1289
 
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3*NPHAS
 
1315
 * Size of DUDXYZ and WDUDXY = n_ghost_cells*3*3
1290
1316
 *----------------------------------------------------------------------------*/
1291
1317
 
1292
1318
void
1293
 
CS_PROCF (peinu2, PEINU2)(const cs_int_t    *iphas,
1294
 
                          cs_real_t         *dudxyz)
 
1319
CS_PROCF (peinu2, PEINU2)(cs_real_t         *dudxyz)
1295
1320
{
1296
1321
  cs_int_t  i, t_id, rank_id, shift;
1297
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
1322
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
1298
1323
  fvm_periodicity_type_t  perio_type;
1299
1324
  cs_real_t  matrix[3][4];
1300
1325
 
1303
1328
 
1304
1329
  const cs_int_t  n_transforms = mesh->n_transforms;
1305
1330
  const cs_int_t  n_ghost_cells = mesh->n_ghost_cells;
1306
 
  const cs_int_t  phase_id = *iphas - 1;
1307
1331
  const cs_int_t  stride = 3 * n_ghost_cells;
1308
1332
  const fvm_periodicity_t  *periodicity = mesh->periodicity;
1309
1333
 
1310
 
  /* Macro pour la position au sein d'un tableau */
1311
 
 
1312
 
#define GET_ID1(i, j, k) ( \
1313
 
   i + n_ghost_cells*j + stride*k + 3*stride*phase_id)
1314
 
 
1315
 
  if (mesh->halo_type == CS_HALO_N_TYPES)
 
1334
  /* Macro for position inside an array */
 
1335
 
 
1336
#define GET_ID1(i, j, k) ( i + n_ghost_cells*j + stride*k )
 
1337
 
 
1338
  if (mesh->halo_type == CS_HALO_N_TYPES || halo == NULL)
1316
1339
    return;
1317
1340
 
1318
 
  assert(halo != NULL);
1319
 
 
1320
1341
  /* Compute the new cell centers through periodicity */
1321
1342
 
1322
1343
  for (t_id = 0; t_id < n_transforms; t_id++) {
1401
1422
 * *****************
1402
1423
 *
1403
1424
 * INTEGER          ISOU          : -> : component of the Reynolds stress tensor
1404
 
 * INTEGER          IPHAS         : -> : current phase number
1405
1425
 * DOUBLE PRECISION DRDXYZ        : -> : gradient of the Reynolds stress tensor
1406
1426
 *                                       for ghost cells and for an explicit
1407
1427
 *                                       treatment of the periodicity.
1408
1428
 * DOUBLE PRECISION W1..3(NCELET) : -  : working buffers
1409
1429
 *
1410
 
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3*NPHAS
 
1430
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3
1411
1431
 *----------------------------------------------------------------------------*/
1412
1432
 
1413
1433
void
1414
1434
CS_PROCF (peinr1, PEINR1)(const cs_int_t    *isou,
1415
 
                          const cs_int_t    *iphas,
1416
1435
                          cs_real_t         *drdxyz,
1417
1436
                          cs_real_t          w1[],
1418
1437
                          cs_real_t          w2[],
1422
1441
 
1423
1442
  const cs_int_t  n_ghost_cells = mesh->n_ghost_cells;
1424
1443
  const cs_int_t  comp_id = *isou - 1;
1425
 
  const cs_int_t  phase_id = *iphas - 1;
1426
1444
  const cs_int_t  strid_v = 2 * n_ghost_cells * 3;
1427
 
  const cs_int_t  strid_p = strid_v * 3 * phase_id;
1428
1445
  const cs_int_t  strid_c = n_ghost_cells * comp_id;
1429
1446
 
1430
 
  _peinur1(strid_c, strid_v, strid_p, drdxyz, w1, w2, w3);
 
1447
  _peinur1(strid_c, strid_v, drdxyz, w1, w2, w3);
1431
1448
}
1432
1449
 
1433
1450
/*----------------------------------------------------------------------------
1438
1455
 * SUBROUTINE PEINR2 (VAR)
1439
1456
 * *****************
1440
1457
 *
1441
 
 * INTEGER          IPHAS         :  -> : current phase number
1442
1458
 * DOUBLE PRECISION DRDXYZ        :  -> : gradient of the Reynolds stress tensor
1443
1459
 *                                       for ghost cells and for an explicit
1444
1460
 *                                       treatment of the periodicity.
1445
1461
 *
1446
 
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3*NPHAS
 
1462
 * Size of DRDXYZ and WDRDXY = n_ghost_cells*6*3
1447
1463
 *----------------------------------------------------------------------------*/
1448
1464
 
1449
1465
void
1450
 
CS_PROCF (peinr2, PEINR2)(const cs_int_t    *iphas,
1451
 
                          cs_real_t         *drdxyz)
 
1466
CS_PROCF (peinr2, PEINR2)(cs_real_t         *drdxyz)
1452
1467
{
1453
1468
  cs_int_t  i, i1, i2, j, k, l, m, rank_id, shift, t_id;
1454
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
1469
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
1455
1470
  fvm_periodicity_type_t  perio_type;
1456
1471
 
1457
1472
  cs_real_t  matrix[3][4];
1465
1480
  const cs_int_t  n_transforms = mesh->n_transforms;
1466
1481
  const cs_int_t  n_ghost_cells = mesh->n_ghost_cells;
1467
1482
  const cs_int_t  stride = 2 * 3 * n_ghost_cells;
1468
 
  const cs_int_t  phase_id = *iphas - 1;
1469
1483
  const cs_halo_type_t  halo_type = mesh->halo_type;
1470
1484
  const fvm_periodicity_t  *periodicity = mesh->periodicity;
1471
1485
 
1472
 
#define GET_ID2(elt_id, i, j)   ( \
1473
 
 elt_id + n_ghost_cells*i + stride*j + 3*stride*phase_id)
 
1486
#define GET_ID2(elt_id, i, j)   ( elt_id + n_ghost_cells*i + stride*j )
1474
1487
 
1475
 
  if (halo_type == CS_HALO_N_TYPES)
 
1488
  if (mesh->halo_type == CS_HALO_N_TYPES || halo == NULL)
1476
1489
    return;
1477
1490
 
1478
1491
  assert(halo != NULL);
1796
1809
                     cs_real_t       *coords)
1797
1810
{
1798
1811
  cs_int_t  i, rank_id, t_id, shift;
1799
 
  cs_int_t  start_std, start_ext, end_std, end_ext, length;
 
1812
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
1800
1813
 
1801
1814
  cs_real_t  matrix[3][4];
1802
1815
 
1876
1889
                       cs_real_t        var[])
1877
1890
{
1878
1891
  cs_int_t  i, rank_id, shift, t_id;
1879
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
1892
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
1880
1893
 
1881
1894
  fvm_periodicity_type_t  perio_type = FVM_PERIODICITY_NULL;
1882
1895
 
1941
1954
}
1942
1955
 
1943
1956
/*----------------------------------------------------------------------------
 
1957
 * Synchronize values for a real vector (interleaved) between periodic cells.
 
1958
 *
 
1959
 * parameters:
 
1960
 *   halo      <-> halo associated with variable to synchronize
 
1961
 *   sync_mode --> kind of halo treatment (standard or extended)
 
1962
 *   var       <-> vector to update
 
1963
 *   incvar    <-- specifies the increment for the elements of var
 
1964
 *----------------------------------------------------------------------------*/
 
1965
 
 
1966
void
 
1967
cs_perio_sync_var_vect(const cs_halo_t  *halo,
 
1968
                       cs_halo_type_t    sync_mode,
 
1969
                       cs_real_t         var[],
 
1970
                       int               incvar)
 
1971
{
 
1972
  cs_int_t  i, rank_id, shift, t_id;
 
1973
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
 
1974
 
 
1975
  cs_real_t matrix[3][4];
 
1976
 
 
1977
  fvm_periodicity_type_t  perio_type = FVM_PERIODICITY_NULL;
 
1978
 
 
1979
  const cs_int_t  n_transforms = halo->n_transforms;
 
1980
  const cs_int_t  n_elts   = halo->n_local_elts;
 
1981
  const fvm_periodicity_t  *periodicity = cs_glob_mesh->periodicity;
 
1982
 
 
1983
  if (sync_mode == CS_HALO_N_TYPES)
 
1984
    return;
 
1985
 
 
1986
  assert(halo != NULL);
 
1987
 
 
1988
  _test_halo_compatibility(halo);
 
1989
 
 
1990
  for (t_id = 0; t_id < n_transforms; t_id++) {
 
1991
 
 
1992
    shift = 4 * halo->n_c_domains * t_id;
 
1993
 
 
1994
    perio_type = fvm_periodicity_get_type(periodicity, t_id);
 
1995
 
 
1996
    if (perio_type >= FVM_PERIODICITY_ROTATION) {
 
1997
 
 
1998
      fvm_periodicity_get_matrix(periodicity, t_id, matrix);
 
1999
 
 
2000
      for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
 
2001
 
 
2002
        start_std = n_elts + halo->perio_lst[shift + 4*rank_id];
 
2003
        length = halo->perio_lst[shift + 4*rank_id + 1];
 
2004
        end_std = start_std + length;
 
2005
 
 
2006
        if (sync_mode == CS_HALO_EXTENDED) {
 
2007
 
 
2008
          start_ext = n_elts + halo->perio_lst[shift + 4*rank_id + 2];
 
2009
          length = halo->perio_lst[shift + 4*rank_id + 3];
 
2010
          end_ext = start_ext + length;
 
2011
 
 
2012
        }
 
2013
 
 
2014
        for (i = start_std; i < end_std; i++)
 
2015
          _apply_vector_rotation_i(matrix, var + i*incvar);
 
2016
 
 
2017
        if (sync_mode == CS_HALO_EXTENDED) {
 
2018
 
 
2019
          for (i = start_ext; i < end_ext; i++)
 
2020
            _apply_vector_rotation_i(matrix, var + i*incvar);
 
2021
 
 
2022
        }
 
2023
 
 
2024
      } /* End of loop on ranks */
 
2025
 
 
2026
    } /* End of the treatment of rotation */
 
2027
 
 
2028
  } /* End of loop on transformations */
 
2029
}
 
2030
 
 
2031
/*----------------------------------------------------------------------------
1944
2032
 * Synchronize values for a real vector between periodic cells.
1945
2033
 *
1946
2034
 * parameters:
1954
2042
 *----------------------------------------------------------------------------*/
1955
2043
 
1956
2044
void
1957
 
cs_perio_sync_var_vect(const cs_halo_t *halo,
1958
 
                       cs_halo_type_t   sync_mode,
1959
 
                       cs_perio_rota_t  rota_mode,
1960
 
                       cs_real_t        var_x[],
1961
 
                       cs_real_t        var_y[],
1962
 
                       cs_real_t        var_z[])
 
2045
cs_perio_sync_var_vect_ni(const cs_halo_t *halo,
 
2046
                          cs_halo_type_t   sync_mode,
 
2047
                          cs_perio_rota_t  rota_mode,
 
2048
                          cs_real_t        var_x[],
 
2049
                          cs_real_t        var_y[],
 
2050
                          cs_real_t        var_z[])
1963
2051
{
1964
2052
  cs_int_t  i, rank_id, shift, t_id;
1965
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
2053
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
1966
2054
  cs_real_t  x_in, y_in, z_in;
1967
2055
 
1968
2056
  cs_real_t matrix[3][4];
2103
2191
 *----------------------------------------------------------------------------*/
2104
2192
 
2105
2193
void
2106
 
cs_perio_sync_var_tens(const cs_halo_t *halo,
2107
 
                       cs_halo_type_t   sync_mode,
2108
 
                       cs_real_t        var11[],
2109
 
                       cs_real_t        var12[],
2110
 
                       cs_real_t        var13[],
2111
 
                       cs_real_t        var21[],
2112
 
                       cs_real_t        var22[],
2113
 
                       cs_real_t        var23[],
2114
 
                       cs_real_t        var31[],
2115
 
                       cs_real_t        var32[],
2116
 
                       cs_real_t        var33[])
 
2194
cs_perio_sync_var_tens_ni(const cs_halo_t *halo,
 
2195
                          cs_halo_type_t   sync_mode,
 
2196
                          cs_real_t        var11[],
 
2197
                          cs_real_t        var12[],
 
2198
                          cs_real_t        var13[],
 
2199
                          cs_real_t        var21[],
 
2200
                          cs_real_t        var22[],
 
2201
                          cs_real_t        var23[],
 
2202
                          cs_real_t        var31[],
 
2203
                          cs_real_t        var32[],
 
2204
                          cs_real_t        var33[])
2117
2205
{
2118
2206
  cs_int_t  i, rank_id, shift, t_id;
2119
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
2207
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
2120
2208
  cs_real_t  v11, v12, v13, v21, v22, v23, v31, v32, v33;
2121
2209
 
2122
2210
  cs_real_t  matrix[3][4];
2226
2314
}
2227
2315
 
2228
2316
/*----------------------------------------------------------------------------
 
2317
 * Synchronize values for a real tensor (interleaved) between periodic cells.
 
2318
 *
 
2319
 * parameters:
 
2320
 *   halo      <-> halo associated with variable to synchronize
 
2321
 *   sync_mode --> kind of halo treatment (standard or extended)
 
2322
 *   var       <-> tensor to update
 
2323
 *----------------------------------------------------------------------------*/
 
2324
 
 
2325
void
 
2326
cs_perio_sync_var_tens(const cs_halo_t *halo,
 
2327
                       cs_halo_type_t   sync_mode,
 
2328
                       cs_real_t        var[])
 
2329
{
 
2330
  cs_int_t  i, rank_id, shift, t_id;
 
2331
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
 
2332
  cs_real_t  v11, v12, v13, v21, v22, v23, v31, v32, v33;
 
2333
 
 
2334
  cs_real_t  matrix[3][4];
 
2335
 
 
2336
  fvm_periodicity_type_t  perio_type = FVM_PERIODICITY_NULL;
 
2337
 
 
2338
  const cs_int_t  n_transforms = halo->n_transforms;
 
2339
  const cs_int_t  n_elts   = halo->n_local_elts;
 
2340
  const fvm_periodicity_t *periodicity = cs_glob_mesh->periodicity;
 
2341
 
 
2342
  if (sync_mode == CS_HALO_N_TYPES)
 
2343
    return;
 
2344
 
 
2345
  assert(halo != NULL);
 
2346
 
 
2347
  _test_halo_compatibility(halo);
 
2348
 
 
2349
  for (t_id = 0; t_id < n_transforms; t_id++) {
 
2350
 
 
2351
    shift = 4 * halo->n_c_domains * t_id;
 
2352
 
 
2353
    perio_type = fvm_periodicity_get_type(periodicity, t_id);
 
2354
 
 
2355
    if (perio_type >= FVM_PERIODICITY_ROTATION) {
 
2356
 
 
2357
      fvm_periodicity_get_matrix(periodicity, t_id, matrix);
 
2358
 
 
2359
      for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
 
2360
 
 
2361
        start_std = halo->perio_lst[shift + 4*rank_id];
 
2362
        length = halo->perio_lst[shift + 4*rank_id + 1];
 
2363
        end_std = start_std + length;
 
2364
 
 
2365
        if (sync_mode == CS_HALO_EXTENDED) {
 
2366
 
 
2367
          start_ext = halo->perio_lst[shift + 4*rank_id + 2];
 
2368
          length = halo->perio_lst[shift + 4*rank_id + 3];
 
2369
          end_ext = start_ext + length;
 
2370
 
 
2371
        }
 
2372
 
 
2373
        for (i = start_std; i < end_std; i++) {
 
2374
 
 
2375
          v11 = var[0 + 3*0 + 9*(n_elts + i)];
 
2376
          v12 = var[0 + 3*1 + 9*(n_elts + i)];
 
2377
          v13 = var[0 + 3*2 + 9*(n_elts + i)];
 
2378
          v21 = var[1 + 3*0 + 9*(n_elts + i)];
 
2379
          v22 = var[1 + 3*1 + 9*(n_elts + i)];
 
2380
          v23 = var[1 + 3*2 + 9*(n_elts + i)];
 
2381
          v31 = var[2 + 3*0 + 9*(n_elts + i)];
 
2382
          v32 = var[2 + 3*1 + 9*(n_elts + i)];
 
2383
          v33 = var[2 + 3*2 + 9*(n_elts + i)];
 
2384
 
 
2385
          _apply_tensor_rotation(matrix,
 
2386
                                 v11, v12, v13, v21, v22, v23,
 
2387
                                 v31, v32, v33,
 
2388
                                 &var[0 + 3*0 + 9*(n_elts + i)],
 
2389
                                 &var[0 + 3*1 + 9*(n_elts + i)],
 
2390
                                 &var[0 + 3*2 + 9*(n_elts + i)],
 
2391
                                 &var[1 + 3*0 + 9*(n_elts + i)],
 
2392
                                 &var[1 + 3*1 + 9*(n_elts + i)],
 
2393
                                 &var[1 + 3*2 + 9*(n_elts + i)],
 
2394
                                 &var[2 + 3*0 + 9*(n_elts + i)],
 
2395
                                 &var[2 + 3*1 + 9*(n_elts + i)],
 
2396
                                 &var[2 + 3*2 + 9*(n_elts + i)]);
 
2397
 
 
2398
        }
 
2399
 
 
2400
        if (sync_mode == CS_HALO_EXTENDED) {
 
2401
 
 
2402
          for (i = start_ext; i < end_ext; i++) {
 
2403
 
 
2404
            v11 = var[0 + 3*0 + 9*(n_elts + i)];
 
2405
            v12 = var[0 + 3*1 + 9*(n_elts + i)];
 
2406
            v13 = var[0 + 3*2 + 9*(n_elts + i)];
 
2407
            v21 = var[1 + 3*0 + 9*(n_elts + i)];
 
2408
            v22 = var[1 + 3*1 + 9*(n_elts + i)];
 
2409
            v23 = var[1 + 3*2 + 9*(n_elts + i)];
 
2410
            v31 = var[2 + 3*0 + 9*(n_elts + i)];
 
2411
            v32 = var[2 + 3*1 + 9*(n_elts + i)];
 
2412
            v33 = var[2 + 3*2 + 9*(n_elts + i)];
 
2413
 
 
2414
            _apply_tensor_rotation(matrix,
 
2415
                                   v11, v12, v13, v21, v22, v23,
 
2416
                                   v31, v32, v33,
 
2417
                                   &var[0 + 3*0 + 9*(n_elts + i)],
 
2418
                                   &var[0 + 3*1 + 9*(n_elts + i)],
 
2419
                                   &var[0 + 3*2 + 9*(n_elts + i)],
 
2420
                                   &var[1 + 3*0 + 9*(n_elts + i)],
 
2421
                                   &var[1 + 3*1 + 9*(n_elts + i)],
 
2422
                                   &var[1 + 3*2 + 9*(n_elts + i)],
 
2423
                                   &var[2 + 3*0 + 9*(n_elts + i)],
 
2424
                                   &var[2 + 3*1 + 9*(n_elts + i)],
 
2425
                                   &var[2 + 3*2 + 9*(n_elts + i)]);
 
2426
 
 
2427
          }
 
2428
 
 
2429
        } /* End of the treatment of rotation */
 
2430
 
 
2431
      } /* End if halo is extended */
 
2432
 
 
2433
    } /* End of loop on ranks */
 
2434
 
 
2435
  } /* End of loop on transformations for the local rank */
 
2436
}
 
2437
 
 
2438
/*----------------------------------------------------------------------------
2229
2439
 * Synchronize values for a real diagonal tensor between periodic cells.
2230
2440
 *
2231
2441
 * We only know the diagonal of the tensor.
2239
2449
 *----------------------------------------------------------------------------*/
2240
2450
 
2241
2451
void
2242
 
cs_perio_sync_var_diag(const cs_halo_t *halo,
2243
 
                       cs_halo_type_t   sync_mode,
2244
 
                       cs_real_t        var11[],
2245
 
                       cs_real_t        var22[],
2246
 
                       cs_real_t        var33[])
 
2452
cs_perio_sync_var_diag_ni(const cs_halo_t *halo,
 
2453
                          cs_halo_type_t   sync_mode,
 
2454
                          cs_real_t        var11[],
 
2455
                          cs_real_t        var22[],
 
2456
                          cs_real_t        var33[])
2247
2457
{
2248
2458
  cs_int_t  i, rank_id, shift, t_id;
2249
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
2459
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
2250
2460
  cs_real_t  v11, v22, v33;
2251
2461
  cs_real_t  matrix[3][4];
2252
2462
 
2331
2541
}
2332
2542
 
2333
2543
/*----------------------------------------------------------------------------
 
2544
 * Synchronize values for a real diagonal tensor (interleaved)
 
2545
 * between periodic cells.
 
2546
 *
 
2547
 * We only know the interleaved diagonal of the tensor.
 
2548
 *
 
2549
 * parameters:
 
2550
 *   halo      <-> halo associated with variable to synchronize
 
2551
 *   sync_mode --> kind of halo treatment (standard or extended)
 
2552
 *   var       <-> diagonal tensor to update
 
2553
 *----------------------------------------------------------------------------*/
 
2554
 
 
2555
void
 
2556
cs_perio_sync_var_diag(const cs_halo_t *halo,
 
2557
                       cs_halo_type_t   sync_mode,
 
2558
                       cs_real_t        var[])
 
2559
{
 
2560
  cs_int_t  i, rank_id, shift, t_id;
 
2561
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
 
2562
  cs_real_t  v11, v22, v33;
 
2563
  cs_real_t  matrix[3][4];
 
2564
 
 
2565
  fvm_periodicity_type_t  perio_type = FVM_PERIODICITY_NULL;
 
2566
 
 
2567
  const cs_int_t  n_transforms = halo->n_transforms;
 
2568
  const cs_int_t  n_elts = halo->n_local_elts;
 
2569
  const fvm_periodicity_t *periodicity = cs_glob_mesh->periodicity;
 
2570
 
 
2571
  if (sync_mode == CS_HALO_N_TYPES)
 
2572
    return;
 
2573
 
 
2574
  assert(halo != NULL);
 
2575
 
 
2576
  _test_halo_compatibility(halo);
 
2577
 
 
2578
  for (t_id = 0; t_id < n_transforms; t_id++) {
 
2579
 
 
2580
    shift = 4 * halo->n_c_domains * t_id;
 
2581
 
 
2582
    perio_type = fvm_periodicity_get_type(periodicity, t_id);
 
2583
 
 
2584
    if (perio_type >= FVM_PERIODICITY_ROTATION) {
 
2585
 
 
2586
      fvm_periodicity_get_matrix(periodicity, t_id, matrix);
 
2587
 
 
2588
      for (rank_id = 0; rank_id < halo->n_c_domains; rank_id++) {
 
2589
 
 
2590
        start_std = halo->perio_lst[shift + 4*rank_id];
 
2591
        length = halo->perio_lst[shift + 4*rank_id + 1];
 
2592
        end_std = start_std + length;
 
2593
 
 
2594
        if (sync_mode == CS_HALO_EXTENDED) {
 
2595
 
 
2596
          start_ext = halo->perio_lst[shift + 4*rank_id + 2];
 
2597
          length = halo->perio_lst[shift + 4*rank_id + 3];
 
2598
          end_ext = start_ext + length;
 
2599
 
 
2600
        }
 
2601
 
 
2602
        for (i = start_std; i < end_std; i++) {
 
2603
 
 
2604
          v11 = var[0 + 3*(n_elts + i)];
 
2605
          v22 = var[1 + 3*(n_elts + i)];
 
2606
          v33 = var[2 + 3*(n_elts + i)];
 
2607
 
 
2608
          _apply_tensor_rotation(matrix,
 
2609
                                 v11,   0,   0,   0, v22,   0,
 
2610
                                 0,   0, v33,
 
2611
                                 &var[0 + 3*(n_elts + i)],
 
2612
                                 NULL,
 
2613
                                 NULL,
 
2614
                                 NULL,
 
2615
                                 &var[1 + 3*(n_elts + i)],
 
2616
                                 NULL,
 
2617
                                 NULL,
 
2618
                                 NULL,
 
2619
                                 &var[2 + 3*(n_elts + i)]);
 
2620
 
 
2621
        }
 
2622
 
 
2623
        if (sync_mode == CS_HALO_EXTENDED) {
 
2624
 
 
2625
          for (i = start_ext; i < end_ext; i++) {
 
2626
 
 
2627
            v11 = var[0 + 3*(n_elts + i)];
 
2628
            v22 = var[1 + 3*(n_elts + i)];
 
2629
            v33 = var[2 + 3*(n_elts + i)];
 
2630
 
 
2631
            _apply_tensor_rotation(matrix,
 
2632
                                   v11,   0,   0,   0, v22,   0,
 
2633
                                   0,   0, v33,
 
2634
                                   &var[0 + 3*(n_elts + i)],
 
2635
                                   NULL,
 
2636
                                   NULL,
 
2637
                                   NULL,
 
2638
                                   &var[1 + 3*(n_elts + i)],
 
2639
                                   NULL,
 
2640
                                   NULL,
 
2641
                                   NULL,
 
2642
                                   &var[2 + 3*(n_elts + i)]);
 
2643
 
 
2644
          }
 
2645
 
 
2646
        } /* End if halo is extended */
 
2647
 
 
2648
      } /* End of loop on ranks */
 
2649
 
 
2650
    } /* End of the treatment of rotation */
 
2651
 
 
2652
  } /* End of loop on transformations */
 
2653
}
 
2654
/*----------------------------------------------------------------------------
2334
2655
 * Update global halo backup buffer size so as to be usable with a given halo.
2335
2656
 *
2336
2657
 * This function should be called at the end of any halo creation,
2414
2735
                            const cs_real_t    var[])
2415
2736
{
2416
2737
  cs_int_t  i, rank_id, shift, t_id;
2417
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
2738
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
2418
2739
 
2419
2740
  cs_mesh_t  *mesh = cs_glob_mesh;
2420
2741
  fvm_periodicity_type_t  perio_type = FVM_PERIODICITY_NULL;
2495
2816
                               cs_real_t          var[])
2496
2817
{
2497
2818
  cs_int_t  i, rank_id, shift, t_id;
2498
 
  cs_int_t  start_std, end_std, length, start_ext, end_ext;
 
2819
  cs_int_t  start_std = 0, end_std = 0, length = 0, start_ext = 0, end_ext = 0;
2499
2820
 
2500
2821
  cs_mesh_t  *mesh = cs_glob_mesh;
2501
2822
  fvm_periodicity_type_t  perio_type = FVM_PERIODICITY_NULL;
2570
2891
  return restore_count;
2571
2892
}
2572
2893
 
2573
 
/*----------------------------------------------------------------------------
2574
 
 * Define parameters for building an interface set structure on the main mesh.
2575
 
 *
2576
 
 * parameters:
2577
 
 *   p_n_periodic_lists   <-> pointer to the number of periodic lists (may
2578
 
 *                            be local)
2579
 
 *   p_periodic_num       <-> pointer to the periodicity number (1 to n)
2580
 
 *                            associated with each periodic list (primary
2581
 
 *                            periodicities only)
2582
 
 *   p_n_periodic_couples <-> pointer to the number of periodic couples
2583
 
 *   p_periodic_couples   <-> pointer to the periodic couple list
2584
 
 *----------------------------------------------------------------------------*/
2585
 
 
2586
 
void
2587
 
cs_perio_define_couples(int         *p_n_periodic_lists,
2588
 
                        int         *p_periodic_num[],
2589
 
                        fvm_lnum_t  *p_n_periodic_couples[],
2590
 
                        fvm_gnum_t  **p_periodic_couples[])
2591
 
{
2592
 
 cs_int_t  i, j, k, n_elts;
2593
 
 cs_int_t  face_rank, perio_id, rank_id, fac_id, fac_num;
2594
 
 cs_int_t  shift_key, shift_key_couple, rshift, sshift;
2595
 
 
2596
 
  cs_int_t  n_max_face_vertices = 0, n_face_vertices = 0;
2597
 
  cs_int_t  *send_count = NULL, *recv_count = NULL;
2598
 
  cs_int_t  *send_shift = NULL, *recv_shift = NULL;
2599
 
  cs_int_t  *send_rank_count = NULL, *recv_rank_count = NULL;
2600
 
  cs_int_t  *send_rank_shift = NULL, *recv_rank_shift = NULL;
2601
 
  cs_int_t  *perio_couple_count = NULL, *perio_couple_shift = NULL;
2602
 
  fvm_gnum_t  *send_buffer = NULL, *recv_buffer = NULL;
2603
 
  fvm_gnum_t  *f1_vertices = NULL, *f2_vertices = NULL;
2604
 
 
2605
 
  int  *periodic_num = NULL;
2606
 
  fvm_lnum_t  *n_periodic_couples = NULL;
2607
 
  fvm_gnum_t  **periodic_couples = NULL;
2608
 
 
2609
 
  const cs_int_t  n_ranks = cs_glob_mesh->n_domains;
2610
 
  const cs_int_t  n_perio = cs_glob_mesh->n_init_perio;
2611
 
  const cs_int_t  n_keys = n_ranks * n_perio;
2612
 
  const cs_int_t  local_rank = (cs_glob_rank_id == -1) ? 0:cs_glob_rank_id;
2613
 
  const cs_int_t  *per_face_idx = cs_glob_mesh_builder->per_face_idx;
2614
 
  const cs_int_t  *per_face_lst = cs_glob_mesh_builder->per_face_lst;
2615
 
  const cs_int_t  *per_rank_lst = cs_glob_mesh_builder->per_rank_lst;
2616
 
  const cs_int_t  *face_vtx_idx = cs_glob_mesh->i_face_vtx_idx;
2617
 
  const cs_int_t  *face_vtx_lst = cs_glob_mesh->i_face_vtx_lst;
2618
 
  const fvm_gnum_t  *g_vtx_num = cs_glob_mesh->global_vtx_num;
2619
 
 
2620
 
  /* Allocate and initialize buffers */
2621
 
 
2622
 
  BFT_MALLOC(send_count, n_keys, cs_int_t);
2623
 
  BFT_MALLOC(recv_count, n_keys, cs_int_t);
2624
 
  BFT_MALLOC(perio_couple_count, n_keys, cs_int_t);
2625
 
  BFT_MALLOC(send_shift, n_keys + 1, cs_int_t);
2626
 
  BFT_MALLOC(recv_shift, n_keys + 1, cs_int_t);
2627
 
  BFT_MALLOC(perio_couple_shift, n_keys + 1, cs_int_t);
2628
 
 
2629
 
  for (i = 0; i < n_keys; i++) {
2630
 
    send_count[i] = 0;
2631
 
    recv_count[i] = 0;
2632
 
    send_shift[i] = 0;
2633
 
    recv_shift[i] = 0;
2634
 
    perio_couple_count[i] = 0;
2635
 
    perio_couple_shift[i] = 0;
2636
 
  }
2637
 
 
2638
 
  send_shift[n_keys] = 0;
2639
 
  recv_shift[n_keys] = 0;
2640
 
  perio_couple_shift[n_keys] = 0;
2641
 
 
2642
 
  /* First step: count number of elements to send */
2643
 
 
2644
 
  for (perio_id = 0; perio_id < n_perio; perio_id++) {
2645
 
 
2646
 
    /*
2647
 
      The buffer per_face_lst is composed by pairs of face.
2648
 
      The first element is the local number of the local face. The sign of the
2649
 
      face indicates the sign of the periodicty.
2650
 
      The second element is the local number of the distant face.
2651
 
      The buffer per_rank_lst owns for each distant face the rank id of the
2652
 
      distant rank associated to the face.
2653
 
    */
2654
 
 
2655
 
    for (i = per_face_idx[perio_id]; i < per_face_idx[perio_id + 1]; i++) {
2656
 
 
2657
 
      fac_num = per_face_lst[2*i];
2658
 
      fac_id = CS_ABS(fac_num) - 1;
2659
 
      n_face_vertices = face_vtx_idx[fac_id+1] - face_vtx_idx[fac_id];
2660
 
      n_max_face_vertices = CS_MAX(n_max_face_vertices, n_face_vertices);
2661
 
 
2662
 
      assert(fac_num != 0);
2663
 
 
2664
 
      if (n_ranks > 1)
2665
 
        face_rank = per_rank_lst[i] - 1;
2666
 
      else
2667
 
        face_rank = 0;
2668
 
 
2669
 
      shift_key = n_perio*face_rank + perio_id;
2670
 
 
2671
 
      if (local_rank == face_rank)
2672
 
        send_count[shift_key] += n_face_vertices;
2673
 
 
2674
 
      else
2675
 
        if (fac_num < 0)
2676
 
          send_count[shift_key] += n_face_vertices;
2677
 
 
2678
 
    } /* End of loop on couples */
2679
 
 
2680
 
  } /* End of loop on periodicities */
2681
 
 
2682
 
  /* Exchange number of elements to exchange */
2683
 
 
2684
 
#if defined(HAVE_MPI)
2685
 
  if (n_ranks > 1)
2686
 
    MPI_Alltoall(send_count, n_perio, CS_MPI_INT,
2687
 
                 recv_count, n_perio, CS_MPI_INT, cs_glob_mpi_comm);
2688
 
#endif
2689
 
 
2690
 
  if (n_ranks == 1)
2691
 
    for (i = 0; i < n_keys; i++)
2692
 
      recv_count[i] = send_count[i];
2693
 
 
2694
 
  /* Build perio_couple_count from recv_count */
2695
 
 
2696
 
  for (rank_id = 0; rank_id < n_ranks; rank_id++)
2697
 
    for (perio_id = 0; perio_id < n_perio; perio_id++)
2698
 
      perio_couple_count[n_ranks*perio_id+rank_id]
2699
 
        = recv_count[n_perio*rank_id+perio_id];
2700
 
 
2701
 
  /* Define index */
2702
 
 
2703
 
  for (i = 0; i < n_keys; i++) {
2704
 
 
2705
 
    recv_shift[i+1] = recv_shift[i] + recv_count[i];
2706
 
    send_shift[i+1] = send_shift[i] + send_count[i];
2707
 
    perio_couple_shift[i+1] = perio_couple_shift[i] + perio_couple_count[i];
2708
 
 
2709
 
  }
2710
 
 
2711
 
  assert(perio_couple_shift[n_keys] == recv_shift[n_keys]);
2712
 
 
2713
 
  /* Allocate and initialize buffers for building periodic interface set */
2714
 
 
2715
 
  BFT_MALLOC(periodic_num, n_perio, int);
2716
 
  BFT_MALLOC(n_periodic_couples, n_perio, fvm_lnum_t);
2717
 
  BFT_MALLOC(periodic_couples, n_perio, fvm_gnum_t *);
2718
 
 
2719
 
  for (perio_id = 0; perio_id < n_perio; perio_id++) {
2720
 
 
2721
 
    periodic_num[perio_id] = perio_id+1;
2722
 
    n_periodic_couples[perio_id] = 0;
2723
 
    periodic_couples[perio_id] = NULL;
2724
 
 
2725
 
    for (rank_id = 0; rank_id < n_ranks; rank_id++)
2726
 
      n_periodic_couples[perio_id] +=
2727
 
        perio_couple_count[n_ranks*perio_id + rank_id];
2728
 
 
2729
 
    BFT_MALLOC(periodic_couples[perio_id],
2730
 
               2*n_periodic_couples[perio_id],
2731
 
               fvm_gnum_t);
2732
 
 
2733
 
  } /* End of loop on periodicities */
2734
 
 
2735
 
  /* Build send buffer and the first part of the periodic couples
2736
 
     in global numbering */
2737
 
 
2738
 
  BFT_MALLOC(send_buffer, send_shift[n_keys], fvm_gnum_t);
2739
 
  BFT_MALLOC(recv_buffer, recv_shift[n_keys], fvm_gnum_t);
2740
 
  BFT_MALLOC(f1_vertices, n_max_face_vertices, fvm_gnum_t);
2741
 
  BFT_MALLOC(f2_vertices, n_max_face_vertices, fvm_gnum_t);
2742
 
 
2743
 
  for (i = 0; i < n_keys; i++) {
2744
 
    send_count[i] = 0;
2745
 
    perio_couple_count[i] = 0;
2746
 
  }
2747
 
 
2748
 
  /* Second step: Fill send_buffer and first part of periodic_couples */
2749
 
 
2750
 
  for (perio_id = 0; perio_id < n_perio; perio_id++) {
2751
 
 
2752
 
    fvm_gnum_t  *_periodic_couples = periodic_couples[perio_id];
2753
 
 
2754
 
    for (i = per_face_idx[perio_id]; i < per_face_idx[perio_id+1]; i++) {
2755
 
 
2756
 
      fac_num = per_face_lst[2*i];
2757
 
      fac_id = CS_ABS(fac_num) - 1;
2758
 
      n_face_vertices = face_vtx_idx[fac_id+1] - face_vtx_idx[fac_id];
2759
 
 
2760
 
      if (n_ranks > 1)
2761
 
        face_rank = per_rank_lst[i] - 1;
2762
 
      else
2763
 
        face_rank = 0;
2764
 
 
2765
 
      shift_key = n_perio*face_rank + perio_id;
2766
 
      shift_key_couple = n_ranks*perio_id + face_rank;
2767
 
 
2768
 
      sshift = send_count[shift_key] + send_shift[shift_key];
2769
 
      rshift =  perio_couple_count[shift_key_couple]
2770
 
              + perio_couple_shift[shift_key_couple]
2771
 
              - perio_couple_shift[n_ranks*perio_id];
2772
 
 
2773
 
      /* Define f1 */
2774
 
 
2775
 
      for (k= 0, j = face_vtx_idx[fac_id]-1;
2776
 
           j < face_vtx_idx[fac_id+1]-1; j++, k++) {
2777
 
        if (g_vtx_num != NULL)
2778
 
          f1_vertices[k] = (fvm_gnum_t)g_vtx_num[face_vtx_lst[j]-1];
2779
 
        else
2780
 
          f1_vertices[k] = (fvm_gnum_t)face_vtx_lst[j];
2781
 
      }
2782
 
 
2783
 
      if (local_rank == face_rank) {
2784
 
 
2785
 
        /* Define f2 */
2786
 
 
2787
 
        fac_id = CS_ABS(per_face_lst[2*i+1]) - 1;
2788
 
 
2789
 
        for (k= 0, j = face_vtx_idx[fac_id]-1;
2790
 
             j < face_vtx_idx[fac_id+1]-1; j++, k++) {
2791
 
          if (g_vtx_num != NULL)
2792
 
            f2_vertices[k] = (fvm_gnum_t)g_vtx_num[face_vtx_lst[j]-1];
2793
 
          else
2794
 
            f2_vertices[k] = (fvm_gnum_t)face_vtx_lst[j];
2795
 
        }
2796
 
 
2797
 
        /* Fill buffers */
2798
 
 
2799
 
        if (fac_num < 0) { /* Send f1 and copy f2 in periodic_couples */
2800
 
 
2801
 
          for (k = 0; k < n_face_vertices; k++)
2802
 
            send_buffer[sshift + k] = f1_vertices[k];
2803
 
 
2804
 
          for (k = 0; k < n_face_vertices; k++)
2805
 
            _periodic_couples[2*(rshift + k)] = f2_vertices[k];
2806
 
 
2807
 
        }
2808
 
        else { /* Send f2 and copy f1 in periodic couples */
2809
 
 
2810
 
          for (k = 0; k < n_face_vertices; k++)
2811
 
            send_buffer[sshift + k] = f2_vertices[k];
2812
 
 
2813
 
          for (k = 0; k < n_face_vertices; k++)
2814
 
            _periodic_couples[2*(rshift + k)] = f1_vertices[k];
2815
 
 
2816
 
        }
2817
 
 
2818
 
        send_count[shift_key] += n_face_vertices;
2819
 
        perio_couple_count[shift_key_couple] += n_face_vertices;
2820
 
 
2821
 
      }
2822
 
      else { /* local_rank != face_rank */
2823
 
 
2824
 
        if (fac_num < 0) { /* Fill send_buffer with f1 */
2825
 
 
2826
 
          for (k = 0; k < n_face_vertices; k++)
2827
 
            send_buffer[sshift + k] = f1_vertices[k];
2828
 
 
2829
 
          send_count[shift_key] += n_face_vertices;
2830
 
 
2831
 
        }
2832
 
        else { /* Fill _periodic_couples with f1 */
2833
 
 
2834
 
          for (k = 0; k < n_face_vertices; k++)
2835
 
            _periodic_couples[2*(rshift + k)] = f1_vertices[k];
2836
 
 
2837
 
          perio_couple_count[shift_key_couple] += n_face_vertices;
2838
 
 
2839
 
        }
2840
 
 
2841
 
      } /* Face_rank != local_rank */
2842
 
 
2843
 
    } /* End of loop on couples */
2844
 
 
2845
 
  } /* End of loop on periodicity */
2846
 
 
2847
 
  /* Exchange buffers */
2848
 
 
2849
 
  if (n_ranks > 1) {
2850
 
 
2851
 
    /* count and shift on ranks for MPI_Alltoallv */
2852
 
 
2853
 
    BFT_MALLOC(send_rank_count, n_ranks, cs_int_t);
2854
 
    BFT_MALLOC(recv_rank_count, n_ranks, cs_int_t);
2855
 
    BFT_MALLOC(send_rank_shift, n_ranks, cs_int_t);
2856
 
    BFT_MALLOC(recv_rank_shift, n_ranks, cs_int_t);
2857
 
 
2858
 
    for (i = 0; i < n_ranks; i++) {
2859
 
 
2860
 
      send_rank_count[i] = 0;
2861
 
      recv_rank_count[i] = 0;
2862
 
 
2863
 
      for (perio_id = 0; perio_id < n_perio; perio_id++) {
2864
 
 
2865
 
        send_rank_count[i] += send_count[n_perio*i + perio_id];
2866
 
        recv_rank_count[i] += recv_count[n_perio*i + perio_id];
2867
 
 
2868
 
      } /* End of loop on periodicities */
2869
 
 
2870
 
    } /* End of loop on ranks */
2871
 
 
2872
 
    send_rank_shift[0] = 0;
2873
 
    recv_rank_shift[0] = 0;
2874
 
 
2875
 
    for (i = 0; i < n_ranks - 1; i++) {
2876
 
 
2877
 
      send_rank_shift[i+1] = send_rank_shift[i] + send_rank_count[i];
2878
 
      recv_rank_shift[i+1] = recv_rank_shift[i] + recv_rank_count[i];
2879
 
 
2880
 
    } /* End of loop on ranks */
2881
 
 
2882
 
#if 0 && defined(DEBUG) && !defined(NDEBUG) /* DEBUG */
2883
 
    for (rank_id = 0; rank_id < n_ranks; rank_id++) {
2884
 
      bft_printf("RANK_ID: %d - send_count: %d - send_shift: %d\n",
2885
 
                 rank_id, send_rank_count[rank_id],
2886
 
                 send_rank_shift[rank_id]);
2887
 
      for (i = 0; i < send_rank_count[rank_id]; i++)
2888
 
        bft_printf("\t%5d | %8d | %12u\n",
2889
 
                   i, i + send_rank_shift[rank_id],
2890
 
                   send_buffer[i + send_rank_shift[rank_id]]);
2891
 
      bft_printf_flush();
2892
 
    }
2893
 
#endif
2894
 
 
2895
 
#if defined(HAVE_MPI)
2896
 
    MPI_Alltoallv(send_buffer, send_rank_count, send_rank_shift, FVM_MPI_GNUM,
2897
 
                  recv_buffer, recv_rank_count, recv_rank_shift, FVM_MPI_GNUM,
2898
 
                  cs_glob_mpi_comm);
2899
 
#endif
2900
 
 
2901
 
#if 0 && defined(DEBUG) && !defined(NDEBUG) /* DEBUG */
2902
 
    for (rank_id = 0; rank_id < n_ranks; rank_id++) {
2903
 
 
2904
 
      bft_printf("RANK_ID: %d - recv_count: %d - recv_shift: %d\n",
2905
 
                 rank_id, recv_rank_count[rank_id],
2906
 
                 recv_rank_shift[rank_id]);
2907
 
      for (i = 0; i < recv_rank_count[rank_id]; i++)
2908
 
        bft_printf("\t%5d | %8d | %12u\n",
2909
 
                   i, i + recv_rank_shift[rank_id],
2910
 
                   recv_buffer[i + recv_rank_shift[rank_id]]);
2911
 
      bft_printf_flush();
2912
 
    }
2913
 
#endif
2914
 
 
2915
 
    BFT_FREE(send_rank_count);
2916
 
    BFT_FREE(recv_rank_count);
2917
 
    BFT_FREE(send_rank_shift);
2918
 
    BFT_FREE(recv_rank_shift);
2919
 
 
2920
 
  } /* End if n_ranks > 1 */
2921
 
 
2922
 
  else {
2923
 
 
2924
 
    assert(n_ranks == 1);
2925
 
    assert(send_shift[n_keys] == recv_shift[n_keys]);
2926
 
 
2927
 
    for (i = 0; i < send_shift[n_keys]; i++)
2928
 
      recv_buffer[i] = send_buffer[i];
2929
 
 
2930
 
  }
2931
 
 
2932
 
  /* Memory management */
2933
 
 
2934
 
  BFT_FREE(send_count);
2935
 
  BFT_FREE(send_shift);
2936
 
  BFT_FREE(send_buffer);
2937
 
  BFT_FREE(f1_vertices);
2938
 
  BFT_FREE(f2_vertices);
2939
 
 
2940
 
  /* Finalize periodic couples definition in the global numbering */
2941
 
 
2942
 
  for (perio_id = 0; perio_id < n_perio; perio_id++) {
2943
 
 
2944
 
    fvm_gnum_t  *_periodic_couples = periodic_couples[perio_id];
2945
 
 
2946
 
    for (rank_id = 0; rank_id < n_ranks; rank_id++) {
2947
 
 
2948
 
      shift_key = n_perio*rank_id + perio_id;
2949
 
      shift_key_couple = n_ranks*perio_id + rank_id;
2950
 
 
2951
 
      assert(recv_count[shift_key] == perio_couple_count[shift_key_couple]);
2952
 
      n_elts = recv_count[shift_key];
2953
 
 
2954
 
      sshift =  perio_couple_shift[shift_key_couple]
2955
 
              - perio_couple_shift[n_ranks*perio_id];
2956
 
      rshift = recv_shift[shift_key];
2957
 
 
2958
 
#if 0
2959
 
      bft_printf("\nPERIO: %d - RANK: %d - N_ELTS: %d\n",
2960
 
                 perio_id, rank_id, n_elts);
2961
 
      bft_printf("shift_key: %d, shift_key_couple: %d\n",
2962
 
                 shift_key, shift_key_couple);
2963
 
      bft_printf("SSHIFT: %d - RSHIFT: %d\n", sshift, rshift);
2964
 
      bft_printf_flush();
2965
 
#endif
2966
 
 
2967
 
      for (j = 0; j < n_elts; j++)
2968
 
        _periodic_couples[2*sshift + 2*j+1] = recv_buffer[rshift + j];
2969
 
 
2970
 
    } /* End of loop on ranks */
2971
 
 
2972
 
  } /* End of loop on periodicities */
2973
 
 
2974
 
  BFT_FREE(recv_count);
2975
 
  BFT_FREE(recv_shift);
2976
 
  BFT_FREE(recv_buffer);
2977
 
  BFT_FREE(perio_couple_count);
2978
 
  BFT_FREE(perio_couple_shift);
2979
 
 
2980
 
  /* Values to return. Assign pointers */
2981
 
 
2982
 
  *p_n_periodic_lists = n_perio;
2983
 
  *p_periodic_num = periodic_num;
2984
 
  *p_n_periodic_couples = n_periodic_couples;
2985
 
  *p_periodic_couples = periodic_couples;
2986
 
 
2987
 
#if 0
2988
 
  for (i = 0; i < n_perio; i++) {
2989
 
    cs_int_t  j;
2990
 
    bft_printf("\n\n  Periodicity number: %d\n", periodic_num[i]);
2991
 
    bft_printf("  Number of couples : %d\n", n_periodic_couples[i]);
2992
 
    for (j = 0; j < n_periodic_couples[i]; j++)
2993
 
      bft_printf("%12d --> %12d\n",
2994
 
                 periodic_couples[i][2*j], periodic_couples[i][2*j + 1]);
2995
 
  }
2996
 
#endif
2997
 
 
2998
 
}
2999
 
 
3000
2894
/*----------------------------------------------------------------------------*/
3001
2895
 
3002
2896
END_C_DECLS