00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024 static char const Id[] = "$Id: parteleser.c 2962 2010-02-04 19:50:44Z scholz $\n\n";
00025 static char const Td[] = "$Today: " __FILE__ " " __DATE__ " " __TIME__ " $\n\n";
00026
00027 #define PREALLOC_DG 25
00028
00029 #include "init.h"
00030 #include "util/util.h"
00031
00032 #ifdef METIS
00033 EXTERN_C_BEGIN
00034 #include "metis.h"
00035 EXTERN_C_END
00036 #endif
00037
00038 int MetisPartition(GridData *gdata,int parts)
00039 {
00040 MagparFunctionLogBegin;
00041
00042 int rank,size;
00043 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00044 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00045
00046 if (rank) {
00047 MagparFunctionLogReturn(0);
00048 }
00049
00050 PetscPrintf(PETSC_COMM_WORLD,
00051 "Meshpartitioning using METIS for %i processors\n",
00052 parts
00053 );
00054
00055 #ifdef METIS
00056
00057
00058
00059
00060
00061
00062
00063
00064 int t_edgecut;
00065 int zero=0,two=2;
00066 METIS_PartMeshNodal(
00067 &gdata->n_ele,
00068 &gdata->n_vert,
00069 gdata->elevert,
00070 &two,
00071 &zero,
00072 &parts,
00073 &t_edgecut,
00074 gdata->elenewproc,
00075 gdata->vertnewproc
00076 );
00077 #else
00078 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"Metis library not linked in, but it is required for volume mesh partitioning.\n");
00079 #endif
00080
00081
00082 for (int i=0; i<gdata->n_vert; i++) {
00083 gdata->vertl2g[i]=i;
00084 }
00085 for (int i=0; i<gdata->n_ele; i++) {
00086 gdata->elel2g[i]=i;
00087 }
00088
00089
00090
00091
00092 ascat(&gdata->vertl2g,gdata->n_vert,1,gdata->vertnewproc,parts,0);
00093 ascat(&gdata->elel2g,gdata->n_ele,1,gdata->elenewproc,parts,0);
00094
00095 MagparFunctionLogReturn(0);
00096 }
00097
00098
00099
00100
00101
00102
00103
00104 int DataPartitionElementsSer(GridData *gdata)
00105 {
00106 MagparFunctionLogBegin;
00107
00108 int rank,size;
00109 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(ierr);
00110 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr);
00111
00112 if (rank) {
00113 MagparFunctionLogReturn(0);
00114 }
00115
00116 PetscTruth flg;
00117 int dometispartition;
00118 ierr = PetscOptionsGetInt(PETSC_NULL,"-metispartition",&dometispartition,&flg);CHKERRQ(ierr);
00119 if (flg!=PETSC_TRUE) {
00120 dometispartition=0;
00121 PetscPrintf(PETSC_COMM_WORLD,
00122 "Option -metispartition not found, using default value %i\n",
00123 dometispartition
00124 );
00125 }
00126
00127 if (dometispartition<0) {
00128 SETERRQ(PETSC_ERR_ARG_OUTOFRANGE,"-metispartition invalid! Exiting!\n");
00129 }
00130 else if (dometispartition==0) {
00131 PetscPrintf(PETSC_COMM_WORLD,"METIS not used for mesh partitioning.\n");
00132 MagparFunctionLogReturn(0);
00133 }
00134
00135 if (dometispartition==1) {
00136 dometispartition=size;
00137 }
00138 if (dometispartition==1) {
00139 PetscPrintf(PETSC_COMM_WORLD,"Mesh partitioning not performed since not necessary on one processor.\n");
00140 MagparFunctionLogReturn(0);
00141 }
00142 if (dometispartition<size) {
00143 PetscPrintf(PETSC_COMM_WORLD,"%i partitions for %i processors: some processors unused!\n",dometispartition,size);
00144 }
00145
00146 ierr = MetisPartition(gdata,dometispartition);CHKERRQ(ierr);
00147 ierr = CheckPartition(gdata,dometispartition);CHKERRQ(ierr);
00148
00149
00150
00151 if (dometispartition>size) {
00152 PetscPrintf(PETSC_COMM_WORLD,
00153 "Remapping %i partitions onto %i processors.\n",
00154 dometispartition,size
00155 );
00156 for (int i=0;i<gdata->n_ele;i++) {
00157 gdata->elenewproc[i]= distint(dometispartition,size,gdata->elenewproc[i]);
00158 }
00159 for (int i=0;i<gdata->n_vert;i++) {
00160 gdata->vertnewproc[i]=distint(dometispartition,size,gdata->vertnewproc[i]);
00161 }
00162 }
00163
00164
00165 IS isvert,isele;
00166 ierr = ISCreateGeneral(PETSC_COMM_SELF,gdata->n_vert,gdata->vertl2g,&isvert);CHKERRQ(ierr);
00167 ierr = ISSetPermutation(isvert);CHKERRQ(ierr);
00168 ierr = ISCreateGeneral(PETSC_COMM_SELF,gdata->n_ele,gdata->elel2g,&isele);CHKERRQ(ierr);
00169 ierr = ISSetPermutation(isele);CHKERRQ(ierr);
00170
00171
00172 ierr = PermuteData(gdata,isvert,isele);CHKERRQ(ierr);
00173
00174 ierr = ISDestroy(isvert);CHKERRQ(ierr);
00175 ierr = ISDestroy(isele);CHKERRQ(ierr);
00176
00177 MagparFunctionLogReturn(0);
00178 }