MOAB: Mesh Oriented datABase
(version 5.4.1)
|
00001 /* 00002 * This commgraph_test primary goal is to setup a communication framework between 00003 * two components, based only on a marker associated to the data being migrated 00004 * between the components 00005 * One use case for us is the physics grid in the atmosphere and spectral dynamics grid in the 00006 * atmosphere 00007 * These could be in theory on different sets of PES, although these 2 models are on the same PES 00008 * 00009 * In our case, spectral element data is hooked using the GLOBAL_DOFS tag of 4x4 ids, associated to 00010 * element 00011 * PhysGrid is coming from an AtmPhys.h5m point cloud grid partitioned in 16 00012 * Spectral mesh is our regular wholeATM_T_01.h5m, which is after one time step 00013 * 00014 * phys grid is very sprinkled in the partition, spectral mesh is more compact; For correct 00015 * identification/correspondence in parallel, it would make sense to use boxes for the spectral mesh 00016 * 00017 * We employ our friends the crystal router, in which we use rendezvous algorithm, to set 00018 * up the communication pattern 00019 * 00020 * input: wholeATM_T.h5m file, on 128 procs, the same file that is used by imoab_coupler test 00021 * input: AtmPhys.h5m file, which contains the physics grid, distributed on 16 processes 00022 * input 2: wholeLND.h5m, which is land distributed on 16 processes too 00023 * 00024 * The communication pattern will be established using a rendezvous method, based on the marker 00025 * (in this case, GLOBAL_ID on vertices on phys grid and GLOBAL_DOFS tag on spectral elements) 00026 * 00027 * in the end, we need to modify tag migrate to move data between these types of components, by 00028 * ID 00029 * 00030 * wholeLnd.h5m has holes in the ID space, that we need to account for; 00031 * In the end, this could be used to send data directly from Dyn atm to land; or to lnd on coupler ? 00032 * 00033 * 00034 */ 00035 00036 #include "moab/Core.hpp" 00037 00038 // MPI includes 00039 #include "moab_mpi.h" 00040 #include "moab/ParallelComm.hpp" 00041 #include "MBParallelConventions.h" 00042 00043 #include "moab/iMOAB.h" 00044 #include "TestUtil.hpp" 00045 #include "moab/CpuTimer.hpp" 00046 #include "moab/ProgOptions.hpp" 00047 #include <iostream> 00048 #include <sstream> 00049 00050 #define CHECKIERR( rc, message ) \ 00051 if( 0 != ( rc ) ) \ 00052 { \ 00053 printf( "%s. ErrorCode = %d.\n", message, rc ); \ 00054 CHECK( 0 ); \ 00055 } 00056 00057 using namespace moab; 00058 00059 // #define VERBOSE 00060 00061 // declare some variables outside main method 00062 // easier to pass them around to the test 00063 int ierr; 00064 int rankInGlobalComm, numProcesses; 00065 MPI_Group jgroup; 00066 std::string atmFilename = TestDir + "unittest/wholeATM_T.h5m"; 00067 // on a regular case, 5 ATM 00068 // cmpatm is for atm on atm pes ! it has the spectral mesh 00069 // cmpphys is for atm on atm phys pes ! it has the point cloud , phys grid 00070 00071 // 00072 int rankInAtmComm = -1; 00073 // it is the spectral mesh unique comp id 00074 int cmpatm = 605; // component ids are unique over all pes, and established in advance; 00075 00076 std::string atmPhysFilename = TestDir + "unittest/AtmPhys.h5m"; 00077 std::string atmPhysOutFilename = "outPhys.h5m"; 00078 std::string atmFilename2 = "wholeATM_new.h5m"; 00079 int rankInPhysComm = -1; 00080 // this will be the physics atm com id; it should be actually 5 00081 int physatm = 5; // component ids are unique over all pes, and established in advance; 00082 00083 // int rankInJointComm = -1; 00084 00085 int nghlay = 0; // number of ghost layers for loading the file 00086 00087 std::vector< int > groupTasks; 00088 int startG1 = 0, startG2 = 0, endG1 = numProcesses - 1, endG2 = numProcesses - 1; 00089 int typeA = 1; // spectral mesh, with GLOBAL_DOFS tags on cells 00090 int typeB = 2; // point cloud mesh, with GLOBAL_ID tag on vertices 00091 00092 std::string readopts( "PARALLEL=READ_PART;PARTITION=PARALLEL_PARTITION;PARALLEL_RESOLVE_SHARED_ENTS" ); 00093 std::string readoptsPC( "PARALLEL=READ_PART;PARTITION=PARALLEL_PARTITION" ); 00094 std::string fileWriteOptions( "PARALLEL=WRITE_PART" ); 00095 std::string tagT( "a2oTbot" ); 00096 std::string tagU( "a2oUbot" ); 00097 std::string tagV( "a2oVbot" ); 00098 std::string separ( ":" ); 00099 std::string tagT1( "a2oTbot_1" ); 00100 std::string tagU1( "a2oUbot_1" ); 00101 std::string tagV1( "a2oVbot_1" ); 00102 std::string tagT2( "a2oTbot_2" ); // just one send 00103 00104 int commgraphtest(); 00105 00106 void testspectral_phys() 00107 { 00108 // no changes 00109 commgraphtest(); 00110 } 00111 00112 void testspectral_lnd() 00113 { 00114 // first model is spectral, second is land 00115 atmPhysFilename = TestDir + "unittest/wholeLnd.h5m"; 00116 atmPhysOutFilename = std::string( "outLnd.h5m" ); 00117 atmFilename2 = std::string( "wholeATM_lnd.h5m" ); 00118 commgraphtest(); 00119 } 00120 00121 void testphysatm_lnd() 00122 { 00123 // use for first file the output "outPhys.h5m" from first test 00124 atmFilename = std::string( "outPhys.h5m" ); 00125 atmPhysFilename = std::string( "outLnd.h5m" ); 00126 atmPhysOutFilename = std::string( "physAtm_lnd.h5m" ); 00127 atmFilename2 = std::string( "physBack_lnd.h5m" ); 00128 tagT = tagT1; 00129 tagU = tagU1; 00130 tagV = tagV1; 00131 tagT1 = std::string( "newT" ); 00132 tagT2 = std::string( "newT2" ); 00133 typeA = 2; 00134 commgraphtest(); 00135 } 00136 int commgraphtest() 00137 { 00138 00139 if( !rankInGlobalComm ) 00140 { 00141 std::cout << " first file: " << atmFilename << "\n on tasks : " << startG1 << ":" << endG1 00142 << "\n second file: " << atmPhysFilename << "\n on tasks : " << startG2 << ":" << endG2 << "\n "; 00143 } 00144 00145 // load files on 2 different communicators, groups 00146 // will create a joint comm for rendezvous 00147 MPI_Group atmPEGroup; 00148 groupTasks.resize( numProcesses, 0 ); 00149 for( int i = startG1; i <= endG1; i++ ) 00150 groupTasks[i - startG1] = i; 00151 00152 ierr = MPI_Group_incl( jgroup, endG1 - startG1 + 1, &groupTasks[0], &atmPEGroup ); 00153 CHECKIERR( ierr, "Cannot create atmPEGroup" ) 00154 00155 groupTasks.clear(); 00156 groupTasks.resize( numProcesses, 0 ); 00157 MPI_Group atmPhysGroup; 00158 for( int i = startG2; i <= endG2; i++ ) 00159 groupTasks[i - startG2] = i; 00160 00161 ierr = MPI_Group_incl( jgroup, endG2 - startG2 + 1, &groupTasks[0], &atmPhysGroup ); 00162 CHECKIERR( ierr, "Cannot create atmPhysGroup" ) 00163 00164 // create 2 communicators, one for each group 00165 int ATM_COMM_TAG = 1; 00166 MPI_Comm atmComm; 00167 // atmComm is for atmosphere app; 00168 ierr = MPI_Comm_create_group( MPI_COMM_WORLD, atmPEGroup, ATM_COMM_TAG, &atmComm ); 00169 CHECKIERR( ierr, "Cannot create atmComm" ) 00170 00171 int PHYS_COMM_TAG = 2; 00172 MPI_Comm physComm; 00173 // physComm is for phys atm app 00174 ierr = MPI_Comm_create_group( MPI_COMM_WORLD, atmPhysGroup, PHYS_COMM_TAG, &physComm ); 00175 CHECKIERR( ierr, "Cannot create physComm" ) 00176 00177 // now, create the joint communicator atm physatm 00178 00179 // 00180 MPI_Group joinAtmPhysAtmGroup; 00181 ierr = MPI_Group_union( atmPEGroup, atmPhysGroup, &joinAtmPhysAtmGroup ); 00182 CHECKIERR( ierr, "Cannot create joint atm - phys atm group" ) 00183 int JOIN_COMM_TAG = 5; 00184 MPI_Comm joinComm; 00185 ierr = MPI_Comm_create_group( MPI_COMM_WORLD, joinAtmPhysAtmGroup, JOIN_COMM_TAG, &joinComm ); 00186 CHECKIERR( ierr, "Cannot create joint atm cou communicator" ) 00187 00188 ierr = iMOAB_Initialize( 0, NULL ); // not really needed anything from argc, argv, yet; maybe we should 00189 CHECKIERR( ierr, "Cannot initialize iMOAB" ) 00190 00191 int cmpAtmAppID = -1; 00192 iMOAB_AppID cmpAtmPID = &cmpAtmAppID; // atm 00193 int physAtmAppID = -1; // -1 means it is not initialized 00194 iMOAB_AppID physAtmPID = &physAtmAppID; // phys atm on phys pes 00195 00196 // load atm mesh 00197 if( atmComm != MPI_COMM_NULL ) 00198 { 00199 MPI_Comm_rank( atmComm, &rankInAtmComm ); 00200 ierr = iMOAB_RegisterApplication( "ATM1", &atmComm, &cmpatm, cmpAtmPID ); 00201 CHECKIERR( ierr, "Cannot register ATM App" ) 00202 00203 // load first model 00204 std::string rdopts = readopts; 00205 if( typeA == 2 ) rdopts = readoptsPC; // point cloud 00206 ierr = iMOAB_LoadMesh( cmpAtmPID, atmFilename.c_str(), rdopts.c_str(), &nghlay ); 00207 CHECKIERR( ierr, "Cannot load ATM mesh" ) 00208 } 00209 00210 // load atm phys mesh 00211 if( physComm != MPI_COMM_NULL ) 00212 { 00213 MPI_Comm_rank( physComm, &rankInPhysComm ); 00214 ierr = iMOAB_RegisterApplication( "PhysATM", &physComm, &physatm, physAtmPID ); 00215 CHECKIERR( ierr, "Cannot register PHYS ATM App" ) 00216 00217 // load phys atm mesh all tests this is PC 00218 ierr = iMOAB_LoadMesh( physAtmPID, atmPhysFilename.c_str(), readoptsPC.c_str(), &nghlay ); 00219 CHECKIERR( ierr, "Cannot load Phys ATM mesh" ) 00220 } 00221 00222 if( MPI_COMM_NULL != joinComm ) 00223 { 00224 ierr = iMOAB_ComputeCommGraph( cmpAtmPID, physAtmPID, &joinComm, &atmPEGroup, &atmPhysGroup, &typeA, &typeB, 00225 &cmpatm, &physatm ); 00226 // it will generate parcomm graph between atm and atmPhys models 00227 // 2 meshes, that are distributed in parallel 00228 CHECKIERR( ierr, "Cannot compute comm graph between the two apps " ) 00229 } 00230 00231 if( atmComm != MPI_COMM_NULL ) 00232 { 00233 // call send tag; 00234 std::string tags = tagT + separ + tagU + separ + tagV; 00235 ierr = iMOAB_SendElementTag( cmpAtmPID, tags.c_str(), &joinComm, &physatm ); 00236 CHECKIERR( ierr, "cannot send tag values" ) 00237 } 00238 00239 if( physComm != MPI_COMM_NULL ) 00240 { 00241 // need to define tag storage 00242 std::string tags1 = tagT1 + separ + tagU1 + separ + tagV1 + separ; 00243 int tagType = DENSE_DOUBLE; 00244 int ndof = 1; 00245 if( typeB == 1 ) ndof = 16; 00246 int tagIndex = 0; 00247 ierr = iMOAB_DefineTagStorage( physAtmPID, tagT1.c_str(), &tagType, &ndof, &tagIndex ); 00248 CHECKIERR( ierr, "failed to define the field tag a2oTbot" ); 00249 00250 ierr = iMOAB_DefineTagStorage( physAtmPID, tagU1.c_str(), &tagType, &ndof, &tagIndex ); 00251 CHECKIERR( ierr, "failed to define the field tag a2oUbot" ); 00252 00253 ierr = iMOAB_DefineTagStorage( physAtmPID, tagV1.c_str(), &tagType, &ndof, &tagIndex ); 00254 CHECKIERR( ierr, "failed to define the field tag a2oVbot" ); 00255 00256 ierr = iMOAB_ReceiveElementTag( physAtmPID, tags1.c_str(), &joinComm, &cmpatm ); 00257 CHECKIERR( ierr, "cannot receive tag values" ) 00258 } 00259 00260 // we can now free the sender buffers 00261 if( atmComm != MPI_COMM_NULL ) 00262 { 00263 ierr = iMOAB_FreeSenderBuffers( cmpAtmPID, &physatm ); 00264 CHECKIERR( ierr, "cannot free buffers" ) 00265 } 00266 00267 if( physComm != MPI_COMM_NULL ) 00268 { 00269 ierr = iMOAB_WriteMesh( physAtmPID, atmPhysOutFilename.c_str(), fileWriteOptions.c_str() ); 00270 } 00271 if( physComm != MPI_COMM_NULL ) 00272 { 00273 // send back first tag only 00274 ierr = iMOAB_SendElementTag( physAtmPID, tagT1.c_str(), &joinComm, &cmpatm ); 00275 CHECKIERR( ierr, "cannot send tag values" ) 00276 } 00277 // receive it in a different tag 00278 if( atmComm != MPI_COMM_NULL ) 00279 { 00280 // need to define tag storage 00281 int tagType = DENSE_DOUBLE; 00282 int ndof = 16; 00283 if( typeA == 2 ) ndof = 1; 00284 int tagIndex = 0; 00285 ierr = iMOAB_DefineTagStorage( cmpAtmPID, tagT2.c_str(), &tagType, &ndof, &tagIndex ); 00286 CHECKIERR( ierr, "failed to define the field tag a2oTbot_2" ); 00287 00288 ierr = iMOAB_ReceiveElementTag( cmpAtmPID, tagT2.c_str(), &joinComm, &physatm ); 00289 CHECKIERR( ierr, "cannot receive tag values a2oTbot_2" ) 00290 } 00291 // now send back one tag , into a different tag, and see if we get the same values back 00292 // we can now free the sender buffers 00293 if( physComm != MPI_COMM_NULL ) 00294 { 00295 ierr = iMOAB_FreeSenderBuffers( physAtmPID, &cmpatm ); 00296 CHECKIERR( ierr, "cannot free buffers " ) 00297 } 00298 if( atmComm != MPI_COMM_NULL ) 00299 { 00300 ierr = iMOAB_WriteMesh( cmpAtmPID, atmFilename2.c_str(), fileWriteOptions.c_str() ); 00301 } 00302 00303 // unregister in reverse order 00304 if( physComm != MPI_COMM_NULL ) 00305 { 00306 ierr = iMOAB_DeregisterApplication( physAtmPID ); 00307 CHECKIERR( ierr, "cannot deregister second app model" ) 00308 } 00309 00310 if( atmComm != MPI_COMM_NULL ) 00311 { 00312 ierr = iMOAB_DeregisterApplication( cmpAtmPID ); 00313 CHECKIERR( ierr, "cannot deregister first app model" ) 00314 } 00315 00316 ierr = iMOAB_Finalize(); 00317 CHECKIERR( ierr, "did not finalize iMOAB" ) 00318 00319 // free atm group and comm 00320 if( MPI_COMM_NULL != atmComm ) MPI_Comm_free( &atmComm ); 00321 MPI_Group_free( &atmPEGroup ); 00322 00323 // free atm phys group and comm 00324 if( MPI_COMM_NULL != physComm ) MPI_Comm_free( &physComm ); 00325 MPI_Group_free( &atmPhysGroup ); 00326 00327 // free atm phys group and comm 00328 if( MPI_COMM_NULL != joinComm ) MPI_Comm_free( &joinComm ); 00329 MPI_Group_free( &joinAtmPhysAtmGroup ); 00330 00331 return 0; 00332 } 00333 int main( int argc, char* argv[] ) 00334 { 00335 00336 MPI_Init( &argc, &argv ); 00337 MPI_Comm_rank( MPI_COMM_WORLD, &rankInGlobalComm ); 00338 MPI_Comm_size( MPI_COMM_WORLD, &numProcesses ); 00339 00340 MPI_Comm_group( MPI_COMM_WORLD, &jgroup ); // all processes in global group 00341 00342 // default: load atm on 2 proc, phys grid on 2 procs, establish comm graph, then migrate 00343 // data from atm pes to phys pes, and back 00344 startG1 = 0, startG2 = 0, endG1 = numProcesses - 1, endG2 = numProcesses - 1; 00345 00346 ProgOptions opts; 00347 opts.addOpt< std::string >( "modelA,m", "first model file ", &atmFilename ); 00348 opts.addOpt< int >( "typeA,t", " type of first model ", &typeA ); 00349 00350 opts.addOpt< std::string >( "modelB,n", "second model file", &atmPhysFilename ); 00351 opts.addOpt< int >( "typeB,v", " type of the second model ", &typeB ); 00352 00353 opts.addOpt< int >( "startAtm,a", "start task for first model layout", &startG1 ); 00354 opts.addOpt< int >( "endAtm,b", "end task for first model layout", &endG1 ); 00355 00356 opts.addOpt< int >( "startPhys,c", "start task for second model layout", &startG2 ); 00357 opts.addOpt< int >( "endPhys,d", "end task for second model layout", &endG2 ); 00358 00359 opts.addOpt< std::string >( "output,o", "output filename", &atmPhysOutFilename ); 00360 00361 opts.addOpt< std::string >( "output,o", "output filename", &atmFilename2 ); 00362 00363 opts.parseCommandLine( argc, argv ); 00364 00365 int num_err = 0; 00366 num_err += RUN_TEST( testspectral_phys ); 00367 00368 // 00369 if( argc == 1 ) 00370 { 00371 num_err += RUN_TEST( testspectral_lnd ); 00372 num_err += RUN_TEST( testphysatm_lnd ); 00373 } 00374 MPI_Group_free( &jgroup ); 00375 00376 MPI_Finalize(); 00377 00378 return num_err; 00379 }