MOAB: Mesh Oriented datABase  (version 5.4.1)
scdseq_timing.cpp
Go to the documentation of this file.
00001 /**
00002  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
00003  * storing and accessing finite element mesh data.
00004  *
00005  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
00006  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
00007  * retains certain rights in this software.
00008  *
00009  * This library is free software; you can redistribute it and/or
00010  * modify it under the terms of the GNU Lesser General Public
00011  * License as published by the Free Software Foundation; either
00012  * version 2.1 of the License, or (at your option) any later version.
00013  *
00014  */
00015 
00016 #include "ScdVertexSeq.hpp"
00017 #include "ScdElementSeq.hpp"
00018 #include "EntitySequenceManager.hpp"
00019 #include "EntitySequence.hpp"
00020 #include "moab/Core.hpp"
00021 #include "moab/ReadUtilIface.hpp"
00022 
00023 #include <iostream>
00024 #include <ctime>
00025 
00026 using namespace moab;
00027 
00028 /*
00029 // some timing/memory measurement includes; memory measurement
00030 // doesn't work on linux, so they're commented out for now to avoid
00031 // platform problems
00032 #include <unistd.h>
00033 #include <termios.h>
00034 #include <sys/ioctl.h>
00035 #include <ctype.h>
00036 #include <sys/resource.h>
00037 #include <sys/types.h>
00038 #include <sys/stat.h>
00039 #include <fcntl.h>
00040 */
00041 
00042 int create_3dtri_3_sequences( Core* gMB, const int intervals, EntityHandle* vstart, EntityHandle* estart );
00043 int create_3dtri_ucd_sequences( Core* gMB, const int intervals, EntityHandle* vstart, EntityHandle* estart );
00044 void print_time();
00045 
00046 int main( int argc, char** argv )
00047 {
00048     int errors = 0;
00049 
00050     // first we need to make a new Core
00051     Core* gMB = new Core();
00052 
00053     // get the intervals
00054     if( argc < 2 )
00055     {
00056         std::cout << "Usage: <scdseq_timing> <#intervals> " << std::endl
00057                   << " where #intervals is the number of intervals on each side of each cube." << std::endl;
00058         return 0;
00059     }
00060 
00061     int intervals;
00062     sscanf( argv[1], "%d", &intervals );
00063 
00064     char do_option[8];
00065 
00066     bool do_scd = true, do_ucd = true;
00067 
00068     if( argc > 2 )
00069     {
00070         sscanf( argv[2], "%s", do_option );
00071         if( do_option[0] == 'u' )
00072             do_scd = false;
00073         else if( do_option[0] == 's' )
00074             do_ucd = false;
00075         else
00076         {
00077             std::cout << "Didn't understand input; doing both scd and ucd." << std::endl;
00078         }
00079     }
00080 
00081     EntityHandle estart[3], vstart[3];
00082     int total_elements = intervals * intervals * intervals;
00083     std::vector< EntityHandle > connect;
00084     ErrorCode result;
00085     clock_t start, stop;
00086     float time;
00087     char inp[1];
00088 
00089     // wait for input to get memory reading
00090     std::cout << "Hit any key and return to continue...";
00091     std::cin >> inp;
00092     std::cout << std::endl;
00093 
00094     if( do_scd )
00095     {
00096 
00097         // create structured mesh
00098         errors = create_3dtri_3_sequences( gMB, intervals, vstart, estart );
00099         if( errors != 0 )
00100         {
00101             std::cout << "Problem creating structured sequences." << std::endl;
00102             return errors;
00103         }
00104 
00105         // get connectivity
00106 
00107         start = clock();
00108         for( int j = 0; j < 3; j++ )
00109         {
00110             for( int i = 0; i < total_elements; i++ )
00111             {
00112                 result = gMB->get_connectivity( estart[j] + i, connect );
00113                 if( MB_SUCCESS != result ) break;
00114                 connect.clear();
00115             }
00116         }
00117         stop = clock();
00118         time = static_cast< float >( stop - start ) / CLOCKS_PER_SEC;
00119 
00120         std::cout << "Time to get connectivity for scd mesh of " << 3 * total_elements << " elements: " << time
00121                   << " seconds." << std::endl;
00122 
00123         print_time();
00124         // wait for input to get memory reading
00125         std::cout << "Hit any key and return to continue...";
00126         std::cin >> inp;
00127         std::cout << std::endl;
00128 
00129         // destroy this mesh
00130         delete gMB;
00131     }
00132 
00133     if( do_ucd )
00134     {
00135 
00136         // now do the same thing, only unstructured
00137         gMB = new Core();
00138 
00139         // create the elements
00140         errors = create_3dtri_ucd_sequences( gMB, intervals, vstart, estart );
00141         if( errors != 0 )
00142         {
00143             std::cout << "Problem creating unstructured sequences." << std::endl;
00144             return errors;
00145         }
00146 
00147         // get connectivity
00148         std::vector< EntityHandle > connect;
00149         start = clock();
00150         for( int j = 0; j < 3; j++ )
00151         {
00152             for( int i = 0; i < total_elements; i++ )
00153             {
00154                 result = gMB->get_connectivity( estart[j] + i, connect );
00155                 if( MB_SUCCESS != result ) break;
00156                 connect.clear();
00157             }
00158         }
00159         stop = clock();
00160         time = static_cast< float >( stop - start ) / CLOCKS_PER_SEC;
00161 
00162         std::cout << "Time to get connectivity for ucd mesh of " << 3 * total_elements << " elements: " << time
00163                   << " seconds." << std::endl;
00164 
00165         print_time();
00166         // wait for input to get memory reading
00167         std::cout << "Hit any key and return to continue...";
00168         std::cin >> inp;
00169         std::cout << std::endl;
00170 
00171         // destroy this mesh
00172         delete gMB;
00173     }
00174 }
00175 
00176 int create_3dtri_3_sequences( Core* gMB, const int intervals, EntityHandle* vstart, EntityHandle* estart )
00177 {
00178     // create 3 brick esequences arranged such that the all share a common (tri-valent) edge;
00179     // orient each region similarly to the 2dtri_3_esequences test problem, swept into 3d in the
00180     // positive k direction.  This direction is divided into intervals intervals
00181     //
00182     // intervals and intervals controls the i and j intervals in region 0, intervals follows from
00183     // that; intervals divides the k axis
00184 
00185     // input is 4 interval settings controlling the 4 degrees of freedom on the interfacesp
00186     int errors = 0;
00187 
00188     ScdVertexSeq* vseq[3];
00189     ScdElementSeq* eseq[3];
00190 
00191     // set vseq parametric spaces directly from intervals-4
00192     // use 0-based parameterization on vseq's just for fun, which means we'll have to transform into
00193     // eseq system
00194     HomCoord vseq0_minmax[2] = { HomCoord( 0, 0, 0 ), HomCoord( intervals, intervals, intervals ) };
00195     HomCoord vseq1_minmax[2] = { HomCoord( 0, 0, 0 ), HomCoord( intervals - 1, intervals, intervals ) };
00196     HomCoord vseq2_minmax[2] = { HomCoord( 0, 0, 0 ), HomCoord( intervals - 1, intervals - 1, intervals ) };
00197 
00198     // get the seq manager from gMB
00199     EntitySequenceManager* seq_mgr = gMB->sequence_manager();
00200 
00201     // create three vertex sequences
00202     EntitySequence* dum_seq;
00203     vseq[0] = vseq[1] = vseq[2] = NULL;
00204 
00205     // first vertex sequence
00206     ErrorCode result =
00207         seq_mgr->create_scd_sequence( vseq0_minmax[0], vseq0_minmax[1], MBVERTEX, 1, vstart[0], dum_seq );
00208     if( NULL != dum_seq ) vseq[0] = dynamic_cast< ScdVertexSeq* >( dum_seq );
00209     assert( MB_FAILURE != result && vstart[0] != 0 && dum_seq != NULL && vseq[0] != NULL );
00210 
00211     // second vertex sequence
00212     result = seq_mgr->create_scd_sequence( vseq1_minmax[0], vseq1_minmax[1], MBVERTEX, 1, vstart[1], dum_seq );
00213     if( NULL != dum_seq ) vseq[1] = dynamic_cast< ScdVertexSeq* >( dum_seq );
00214     assert( MB_FAILURE != result && vstart[1] != 0 && dum_seq != NULL && vseq[1] != NULL );
00215 
00216     // third vertex sequence
00217     result = seq_mgr->create_scd_sequence( vseq2_minmax[0], vseq2_minmax[1], MBVERTEX, 1, vstart[2], dum_seq );
00218     if( NULL != dum_seq ) vseq[2] = dynamic_cast< ScdVertexSeq* >( dum_seq );
00219     assert( MB_FAILURE != result && vstart[2] != 0 && dum_seq != NULL && vseq[2] != NULL );
00220 
00221     // now create the three element sequences
00222 
00223     // set eseq parametric spaces directly from intervals-4
00224     // use 0-based parameterization on eseq's just for fun, which means we'll have to transform into
00225     // eseq system
00226     HomCoord eseq0_minmax[2] = { HomCoord( 0, 0, 0 ), HomCoord( intervals, intervals, intervals ) };
00227     HomCoord eseq1_minmax[2] = { HomCoord( 0, 0, 0 ), HomCoord( intervals, intervals, intervals ) };
00228     HomCoord eseq2_minmax[2] = { HomCoord( 0, 0, 0 ), HomCoord( intervals, intervals, intervals ) };
00229 
00230     eseq[0] = eseq[1] = eseq[2] = NULL;
00231 
00232     // create the first element sequence
00233     result = seq_mgr->create_scd_sequence( eseq0_minmax[0], eseq0_minmax[1], MBHEX, 1, estart[0], dum_seq );
00234     if( NULL != dum_seq ) eseq[0] = dynamic_cast< ScdElementSeq* >( dum_seq );
00235     assert( MB_FAILURE != result && estart[0] != 0 && dum_seq != NULL && eseq[0] != NULL );
00236 
00237     // only need to add one vseq to this, unity transform
00238     result = eseq[0]->add_vsequence( vseq[0],
00239                                      // trick: if I know it's going to be unity, just input
00240                                      // 3 sets of equivalent points
00241                                      vseq0_minmax[0], vseq0_minmax[0], vseq0_minmax[0], vseq0_minmax[0],
00242                                      vseq0_minmax[0], vseq0_minmax[0] );
00243 
00244     if( MB_SUCCESS != result )
00245     {
00246         std::cout << "Couldn't add first vsequence to first element sequence in tri-composite 3d eseq." << std::endl;
00247         errors++;
00248     }
00249 
00250     // create the second element sequence
00251     result = seq_mgr->create_scd_sequence( eseq1_minmax[0], eseq1_minmax[1], MBHEX, 1, estart[1], dum_seq );
00252     if( NULL != dum_seq ) eseq[1] = dynamic_cast< ScdElementSeq* >( dum_seq );
00253     assert( MB_FAILURE != result && estart[1] != 0 && dum_seq != NULL && eseq[1] != NULL );
00254 
00255     // add shared side from first vseq to this eseq, with bb to get just the face
00256     result = eseq[1]->add_vsequence( vseq[0],
00257                                      // p1: origin in both systems
00258                                      vseq0_minmax[0], eseq0_minmax[0],
00259                                      // p2: one unit along the shared line (i in one, j in other)
00260                                      vseq0_minmax[0] + HomCoord::unitv[0], eseq0_minmax[0] + HomCoord::unitv[1],
00261                                      // p3: arbitrary
00262                                      vseq0_minmax[0], eseq0_minmax[0],
00263                                      // set bb such that it's the jmin side of vseq
00264                                      true, eseq[1]->min_params(),
00265                                      HomCoord( eseq[1]->min_params().i(), eseq[1]->max_params().j(),
00266                                                eseq[1]->max_params().k() ) );
00267     if( MB_SUCCESS != result )
00268     {
00269         std::cout << "Couldn't add shared vsequence to second element sequence in tri-composite 3d eseq." << std::endl;
00270         errors++;
00271     }
00272 
00273     // add second vseq to this eseq, with different orientation but all of it (no bb input)
00274     result =
00275         eseq[1]->add_vsequence( vseq[1],
00276                                 // p1: origin/i+1 (vseq/eseq)
00277                                 vseq1_minmax[0], eseq1_minmax[0] + HomCoord::unitv[0],
00278                                 // p2: j+1 from p1
00279                                 vseq1_minmax[0] + HomCoord::unitv[1],
00280                                 eseq1_minmax[0] + HomCoord::unitv[0] + HomCoord::unitv[1],
00281                                 // p3: i+1 from p1
00282                                 vseq1_minmax[0] + HomCoord::unitv[0], eseq[1]->min_params() + HomCoord::unitv[0] * 2 );
00283     if( MB_SUCCESS != result )
00284     {
00285         std::cout << "Couldn't add second vseq to second element sequence in tri-composite 3d eseq." << std::endl;
00286         errors++;
00287     }
00288 
00289     // create the third element sequence
00290     result = seq_mgr->create_scd_sequence( eseq2_minmax[0], eseq2_minmax[1], MBHEX, 1, estart[2], dum_seq );
00291     if( NULL != dum_seq ) eseq[2] = dynamic_cast< ScdElementSeq* >( dum_seq );
00292     assert( MB_FAILURE != result && estart[2] != 0 && dum_seq != NULL && eseq[2] != NULL );
00293 
00294     // add shared side from second vseq to this eseq
00295     result = eseq[2]->add_vsequence(
00296         vseq[1],
00297         // p1: origin/j+1 (vseq/eseq)
00298         vseq1_minmax[0], eseq[2]->min_params() + HomCoord::unitv[1],
00299         // p2: i+1/j+2 (vseq/eseq)
00300         vseq1_minmax[0] + HomCoord::unitv[0], eseq[2]->min_params() + HomCoord::unitv[1] * 2,
00301         // p3: arbitrary
00302         vseq1_minmax[0], eseq[2]->min_params() + HomCoord::unitv[1],
00303         // bb input such that we only get one side of eseq parameter space
00304         true, eseq[2]->min_params() + HomCoord::unitv[1],
00305         HomCoord( eseq[2]->min_params().i(), eseq[2]->max_params().j(), eseq[2]->max_params().k() ) );
00306     if( MB_SUCCESS != result )
00307     {
00308         std::cout << "Couldn't add shared vsequence to third element sequence in tri-composite 3d eseq." << std::endl;
00309         errors++;
00310     }
00311 
00312     // add shared side from first vseq to this eseq
00313     result = eseq[2]->add_vsequence( vseq[0],
00314                                      // p1: origin/origin
00315                                      vseq1_minmax[0], eseq2_minmax[0],
00316                                      // p2: j+1/i+1
00317                                      vseq1_minmax[0] + HomCoord::unitv[1], eseq2_minmax[0] + HomCoord::unitv[0],
00318                                      // p3: arbitrary
00319                                      vseq1_minmax[0], eseq2_minmax[0],
00320                                      // bb input such that we only get one side of eseq parameter space
00321                                      true, eseq2_minmax[0],
00322                                      HomCoord( eseq2_minmax[1].i(), eseq2_minmax[0].j(), eseq2_minmax[1].k() ) );
00323     if( MB_SUCCESS != result )
00324     {
00325         std::cout << "Couldn't add left shared vsequence to third element sequence in "
00326                      "tri-composite 3d eseq."
00327                   << std::endl;
00328         errors++;
00329     }
00330 
00331     // add third vseq to this eseq
00332     result = eseq[2]->add_vsequence( vseq[2],
00333                                      // p1: origin/i+1,j+1
00334                                      vseq2_minmax[0], eseq[2]->min_params() + HomCoord::unitv[0] + HomCoord::unitv[1],
00335                                      // p2: i+1 from p1
00336                                      vseq2_minmax[0] + HomCoord::unitv[0],
00337                                      eseq[2]->min_params() + HomCoord::unitv[0] * 2 + HomCoord::unitv[1],
00338                                      // p3: j+1 from p1
00339                                      vseq2_minmax[0] + HomCoord::unitv[1],
00340                                      eseq[2]->min_params() + HomCoord::unitv[0] + HomCoord::unitv[1] * 2 );
00341     if( MB_SUCCESS != result )
00342     {
00343         std::cout << "Couldn't add third vseq to third element sequence in tri-composite 3d eseq." << std::endl;
00344         errors++;
00345     }
00346 
00347     return errors;
00348 }
00349 
00350 int create_3dtri_ucd_sequences( Core* gMB, const int intervals, EntityHandle* vstart, EntityHandle* estart )
00351 {
00352 
00353     ReadUtilIface* readMeshIface;
00354     gMB->query_interface( readMeshIface );
00355 
00356     int num_elements = intervals * intervals * intervals;
00357     int num_verts    = ( intervals + 1 ) * ( intervals + 1 ) * ( intervals + 1 );
00358 
00359     std::vector< double* > arrays;
00360     for( int i = 0; i < 3; i++ )
00361     {
00362         readMeshIface->get_node_coords( 3, num_verts, MB_START_ID, vstart[i], arrays );
00363         arrays.clear();
00364     }
00365 
00366     EntityHandle* conn[3];
00367 
00368     // allocate 3 arrays to initialize connectivity data
00369     for( int i = 0; i < 3; i++ )
00370     {
00371         readMeshIface->get_element_connect( num_elements, 8, MBHEX, 1, estart[i], conn[i] );
00372 
00373         // now, initialize connectivity data to what it should be; just fudge for now
00374         for( int j = 0; j < num_elements * 8; j++ )
00375             conn[i][j] = vstart[i];
00376     }
00377 
00378     return 0;
00379 }
00380 
00381 void print_time()
00382 {
00383     /*
00384       struct rusage r_usage;
00385       float utime, stime;
00386       getrusage(RUSAGE_SELF, &r_usage);
00387 
00388       if (r_usage.ru_maxrss == 0) {
00389           // this machine doesn't return rss - try going to /proc
00390           // print the file name to open
00391         char file_str[4096], dum_str[4096];
00392         int file_ptr = -1, file_len;
00393         file_ptr = open("/proc/self/stat", O_RDONLY);
00394         file_len = read(file_ptr, file_str, sizeof(file_str)-1);
00395         if (file_len == 0) return;
00396         close(file_ptr);
00397         file_str[file_len] = '\0';
00398           // read the preceeding fields and the ones we really want...
00399         int dum_int;
00400         unsigned int dum_uint, vm_size, rss;
00401         static int page_size = getpagesize();
00402         int num_fields = sscanf(file_str,
00403                                 "%d " // pid
00404                                 "%s " // comm
00405                                 "%c " // state
00406                                 "%d %d %d %d %d " // ppid, pgrp, session, tty, tpgid
00407                                 "%u %u %u %u %u " // flags, minflt, cminflt, majflt, cmajflt
00408                                 "%d %d %d %d %d %d " // utime, stime, cutime, cstime, counter,
00409       priority
00410                                 "%u %u " // timeout, itrealvalue
00411                                 "%d " // starttime
00412                                 "%u %u", // vsize, rss
00413                                 &dum_int,
00414                                 dum_str,
00415                                 dum_str,
00416                                 &dum_int, &dum_int, &dum_int, &dum_int, &dum_int,
00417                                 &dum_uint, &dum_uint, &dum_uint, &dum_uint, &dum_uint,
00418                                 &dum_int, &dum_int, &dum_int, &dum_int, &dum_int, &dum_int,
00419                                 &dum_uint, &dum_uint,
00420                                 &dum_int,
00421                                 &vm_size, &rss);
00422         if (num_fields == 24) {
00423           r_usage.ru_maxrss = rss/page_size;
00424           r_usage.ru_idrss = vm_size/page_size;
00425         }
00426       }
00427       utime = (float)r_usage.ru_utime.tv_sec +
00428         ((float)r_usage.ru_utime.tv_usec/1.e6);
00429       stime = (float)r_usage.ru_stime.tv_sec +
00430         ((float)r_usage.ru_stime.tv_usec/1.e6);
00431       static int pagesize = getpagesize();
00432 
00433       std::cout << "Execution time = " << utime+stime << ", max RSS = "
00434                 << r_usage.ru_maxrss*pagesize << " bytes." << std::endl;
00435 
00436     */
00437 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines