EpetraExt Package Browser (Single Doxygen Collection) Development
rect2DMeshGenerator.cpp
Go to the documentation of this file.
00001 /*
00002 //@HEADER
00003 // ***********************************************************************
00004 //
00005 //     EpetraExt: Epetra Extended - Linear Algebra Services Package
00006 //                 Copyright (2011) Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00039 //
00040 // ***********************************************************************
00041 //@HEADER
00042 */
00043 
00044 #include "rect2DMeshGenerator.hpp"
00045 #include "Epetra_IntSerialDenseVector.h"
00046 #include "Epetra_SerialDenseMatrix.h"
00047 #include "Epetra_IntSerialDenseMatrix.h"
00048 #include "Teuchos_Assert.hpp"
00049 #include "Teuchos_FancyOStream.hpp"
00050 #include "Teuchos_RefCountPtr.hpp"
00051 
00052 void GLpApp::rect2DMeshGenerator(
00053   const int                      numProc
00054   ,const int                     procRank
00055   ,const double                  len_x
00056   ,const double                  len_y
00057   ,const int                     local_nx
00058   ,const int                     local_ny
00059   ,const int                     bndy_marker
00060   ,Epetra_IntSerialDenseVector   *ipindx_out
00061   ,Epetra_SerialDenseMatrix      *ipcoords_out
00062   ,Epetra_IntSerialDenseVector   *pindx_out
00063   ,Epetra_SerialDenseMatrix      *pcoords_out
00064   ,Epetra_IntSerialDenseMatrix   *t_out
00065   ,Epetra_IntSerialDenseMatrix   *e_out
00066   ,std::ostream                  *out_arg
00067   ,const bool                    dumpAll
00068   )
00069 {
00070   Teuchos::RefCountPtr<Teuchos::FancyOStream>
00071     out = Teuchos::getFancyOStream(Teuchos::rcp(out_arg,false));
00072   Teuchos::OSTab tab(out);
00073   if(out.get())
00074     *out << "\nGenerating rectangular mesh with len_x="<<len_x<<", len_y="<<len_y<<", local_nx="<<local_nx<<", and local_ny="<<local_ny<<" ...\n";
00075   //
00076   // Validate input
00077   //
00078 #ifdef TEUCHOS_DEBUG
00079   TEUCHOS_TEST_FOR_EXCEPT(len_x <= 0.0);
00080   TEUCHOS_TEST_FOR_EXCEPT(len_y <= 0.0);
00081   TEUCHOS_TEST_FOR_EXCEPT(local_nx <= 0);
00082   TEUCHOS_TEST_FOR_EXCEPT(local_ny <= 0);
00083   TEUCHOS_TEST_FOR_EXCEPT(ipindx_out==NULL);
00084   TEUCHOS_TEST_FOR_EXCEPT(ipcoords_out==NULL);
00085   TEUCHOS_TEST_FOR_EXCEPT(pindx_out==NULL);
00086   TEUCHOS_TEST_FOR_EXCEPT(pcoords_out==NULL);
00087   TEUCHOS_TEST_FOR_EXCEPT(t_out==NULL);
00088   TEUCHOS_TEST_FOR_EXCEPT(e_out==NULL);
00089 #endif
00090   //
00091   // Get local references
00092   //
00093   Epetra_IntSerialDenseVector   &ipindx = *ipindx_out;
00094   Epetra_SerialDenseMatrix      &ipcoords = *ipcoords_out;
00095   Epetra_IntSerialDenseVector   &pindx = *pindx_out;
00096   Epetra_SerialDenseMatrix      &pcoords = *pcoords_out;
00097   Epetra_IntSerialDenseMatrix   &t = *t_out;
00098   Epetra_IntSerialDenseMatrix   &e = *e_out;
00099   //
00100   // Get dimensions
00101   //
00102   const int
00103     numip = ( local_nx + 1 ) * ( local_ny + ( procRank == numProc-1 ? 1 : 0 ) ),
00104     numcp = ( procRank == numProc-1 ? 0 : (local_nx+1) ),
00105     nump  = numip + numcp,
00106     numelems = 2*local_nx*local_ny,
00107     numedges = 2*local_ny + ( procRank==0 ? local_nx : 0 ) + ( procRank==numProc-1 ? local_nx : 0 );
00108   const double
00109     delta_x = len_x / local_nx,
00110     delta_y = (len_y/numProc) / local_ny;
00111   if(out.get()) {
00112     *out
00113       << "\nnumip = " << numip
00114       << "\nnumcp = " << numcp
00115       << "\nnump = " << nump
00116       << "\nnumelems = " << numelems
00117       << "\nnumedges = " << numedges
00118       << "\ndelta_x = " << delta_x
00119       << "\ndelta_y = " << delta_y
00120       << "\n";
00121   }
00122   //
00123   // Resize mesh storage
00124   //
00125   ipindx.Size(numip);
00126   ipcoords.Shape(numip,2);
00127   pindx.Size(nump);
00128   pcoords.Shape(nump,2);
00129   t.Shape(numelems,3);
00130   e.Shape(numedges,3);
00131   //
00132   // Get the global offsets for the local nodes IDs and y coordinates
00133   //
00134   const int     localNodeOffset   = procRank*(local_nx+1)*local_ny;
00135   const double  localYCoordOffset = procRank*delta_y*local_ny;
00136   //
00137   // Set the local node IDs and their cooridinates
00138   //
00139   if(out.get()) *out << "\nGenerating the node list ...\n";
00140   tab.incrTab();
00141   // Owned and shared
00142   int node_i  = 0;
00143   int global_node_id = localNodeOffset+1;
00144   for( int i_y = 0; i_y < local_ny + 1; ++i_y ) {
00145     for( int i_x = 0; i_x < local_nx + 1; ++i_x ) {
00146       pindx(node_i) = global_node_id;
00147       pcoords(node_i,0) = i_x * delta_x;
00148       pcoords(node_i,1) = i_y * delta_y + localYCoordOffset;
00149       if(out.get() && dumpAll)
00150         *out << "node_i="<<node_i<<",global_node_id="<<global_node_id<<",("<<pcoords(node_i,0)<<","<<pcoords(node_i,1)<<")\n";
00151       ++node_i;
00152       ++global_node_id;
00153     }
00154   }
00155   tab.incrTab(-1);
00156   TEUCHOS_TEST_FOR_EXCEPT(node_i != nump);
00157   // Locally owned only
00158   for( int i = 0; i < numip; ++i ) {
00159     ipindx(i) = pindx(i);
00160     ipcoords(i,0) = pcoords(i,0);
00161     ipcoords(i,1) = pcoords(i,1);
00162   }
00163   //
00164   // Set the elements
00165   //
00166   if(out.get()) *out << "\nGenerating the element list ...\n";
00167   tab.incrTab();
00168   int ele_k = 0;
00169   global_node_id = localNodeOffset+1;
00170   for( int i_y = 0; i_y < local_ny; ++i_y ) {
00171     for( int i_x = 0; i_x < local_nx; ++i_x ) {
00172       t(ele_k,0) = global_node_id;
00173       t(ele_k,1) = global_node_id + (local_nx + 1);
00174       t(ele_k,2) = global_node_id + 1;
00175       if(out.get() && dumpAll)
00176         *out << "ele_k="<<ele_k<<",("<<t(ele_k,0)<<","<<t(ele_k,1)<<","<<t(ele_k,2)<<")\n";
00177       ++ele_k;
00178       t(ele_k,0) = global_node_id + 1;
00179       t(ele_k,1) = global_node_id + (local_nx + 1);
00180       t(ele_k,2) = global_node_id + (local_nx + 1) + 1;
00181       if(out.get() && dumpAll)
00182         *out << "ele_k="<<ele_k<<",("<<t(ele_k,0)<<","<<t(ele_k,1)<<","<<t(ele_k,2)<<")\n";
00183       ++ele_k;
00184       ++global_node_id;
00185     }
00186     ++global_node_id;
00187   }
00188   tab.incrTab(-1);
00189   TEUCHOS_TEST_FOR_EXCEPT(ele_k != numelems);
00190   //
00191   // Set the edges
00192   //
00193   int edge_j = 0;
00194   if(procRank==0) {
00195     // Bottom edges
00196     if(out.get()) *out << "\nGenerating the bottom edges ...\n";
00197     tab.incrTab();
00198     global_node_id = localNodeOffset+1;
00199     for( int i_x = 0; i_x < local_nx; ++i_x ) {
00200       e(edge_j,0) = global_node_id;
00201       e(edge_j,1) = global_node_id + 1;
00202       e(edge_j,2) = bndy_marker;
00203       if(out.get() && dumpAll)
00204         *out << "edge_j="<<edge_j<<",("<<e(edge_j,0)<<","<<e(edge_j,1)<<"),"<<e(edge_j,2)<<"\n";
00205       ++edge_j;
00206       global_node_id += 1;
00207     }
00208     tab.incrTab(-1);
00209   }
00210   // Left edges
00211   if(out.get()) *out << "\nGenerating the left edges ...\n";
00212   tab.incrTab();
00213   global_node_id = localNodeOffset+1;
00214   for( int i_y = 0; i_y < local_ny; ++i_y ) {
00215     e(edge_j,0) = global_node_id;
00216     e(edge_j,1) = global_node_id + (local_nx + 1);
00217     e(edge_j,2) = bndy_marker;
00218     if(out.get() && dumpAll)
00219       *out << "edge_j="<<edge_j<<",("<<e(edge_j,0)<<","<<e(edge_j,1)<<"),"<<e(edge_j,2)<<"\n";
00220     ++edge_j;
00221     global_node_id += (local_nx + 1);
00222   }
00223   tab.incrTab(-1);
00224   // Right edges
00225   if(out.get()) *out << "\nGenerating the right edges ...\n";
00226   tab.incrTab();
00227   global_node_id = localNodeOffset + 1 + local_nx;
00228   for( int i_y = 0; i_y < local_ny; ++i_y ) {
00229     e(edge_j,0) = global_node_id;
00230     e(edge_j,1) = global_node_id + (local_nx + 1);
00231     e(edge_j,2) = bndy_marker;
00232     if(out.get() && dumpAll)
00233       *out << "edge_j="<<edge_j<<",("<<e(edge_j,0)<<","<<e(edge_j,1)<<"),"<<e(edge_j,2)<<"\n";
00234     ++edge_j;
00235     global_node_id += (local_nx + 1);
00236   }
00237   tab.incrTab(-1);
00238   if(procRank==numProc-1) {
00239     // Top edges
00240     if(out.get()) *out << "\nGenerating the top edges ...\n";
00241     tab.incrTab();
00242     global_node_id = localNodeOffset+1+(local_nx+1)*local_ny;
00243     for( int i_x = 0; i_x < local_nx; ++i_x ) {
00244       e(edge_j,0) = global_node_id;
00245       e(edge_j,1) = global_node_id + 1;
00246       e(edge_j,2) = bndy_marker;
00247       if(out.get() && dumpAll)
00248         *out << "edge_j="<<edge_j<<",("<<e(edge_j,0)<<","<<e(edge_j,1)<<"),"<<e(edge_j,2)<<"\n";
00249       ++edge_j;
00250       global_node_id += 1;
00251     }
00252     tab.incrTab(-1);
00253   }
00254   TEUCHOS_TEST_FOR_EXCEPT(edge_j != numedges);
00255 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines