45 #include "Epetra_IntSerialDenseVector.h" 46 #include "Epetra_SerialDenseMatrix.h" 47 #include "Epetra_IntSerialDenseMatrix.h" 59 ,
const int bndy_marker
60 ,Epetra_IntSerialDenseVector *ipindx_out
61 ,Epetra_SerialDenseMatrix *ipcoords_out
62 ,Epetra_IntSerialDenseVector *pindx_out
63 ,Epetra_SerialDenseMatrix *pcoords_out
64 ,Epetra_IntSerialDenseMatrix *t_out
65 ,Epetra_IntSerialDenseMatrix *e_out
66 ,std::ostream *out_arg
71 out = Teuchos::getFancyOStream(
Teuchos::rcp(out_arg,
false));
74 *out <<
"\nGenerating rectangular mesh with len_x="<<len_x<<
", len_y="<<len_y<<
", local_nx="<<local_nx<<
", and local_ny="<<local_ny<<
" ...\n";
93 Epetra_IntSerialDenseVector &ipindx = *ipindx_out;
94 Epetra_SerialDenseMatrix &ipcoords = *ipcoords_out;
95 Epetra_IntSerialDenseVector &pindx = *pindx_out;
96 Epetra_SerialDenseMatrix &pcoords = *pcoords_out;
97 Epetra_IntSerialDenseMatrix &t = *t_out;
98 Epetra_IntSerialDenseMatrix &e = *e_out;
103 numip = ( local_nx + 1 ) * ( local_ny + ( procRank == numProc-1 ? 1 : 0 ) ),
104 numcp = ( procRank == numProc-1 ? 0 : (local_nx+1) ),
105 nump = numip + numcp,
106 numelems = 2*local_nx*local_ny,
107 numedges = 2*local_ny + ( procRank==0 ? local_nx : 0 ) + ( procRank==numProc-1 ? local_nx : 0 );
109 delta_x = len_x / local_nx,
110 delta_y = (len_y/numProc) / local_ny;
113 <<
"\nnumip = " << numip
114 <<
"\nnumcp = " << numcp
115 <<
"\nnump = " << nump
116 <<
"\nnumelems = " << numelems
117 <<
"\nnumedges = " << numedges
118 <<
"\ndelta_x = " << delta_x
119 <<
"\ndelta_y = " << delta_y
126 ipcoords.Shape(numip,2);
128 pcoords.Shape(nump,2);
134 const int localNodeOffset = procRank*(local_nx+1)*local_ny;
135 const double localYCoordOffset = procRank*delta_y*local_ny;
139 if(out.
get()) *out <<
"\nGenerating the node list ...\n";
143 int global_node_id = localNodeOffset+1;
144 for(
int i_y = 0; i_y < local_ny + 1; ++i_y ) {
145 for(
int i_x = 0; i_x < local_nx + 1; ++i_x ) {
146 pindx(node_i) = global_node_id;
147 pcoords(node_i,0) = i_x * delta_x;
148 pcoords(node_i,1) = i_y * delta_y + localYCoordOffset;
149 if(out.
get() && dumpAll)
150 *out <<
"node_i="<<node_i<<
",global_node_id="<<global_node_id<<
",("<<pcoords(node_i,0)<<
","<<pcoords(node_i,1)<<
")\n";
158 for(
int i = 0; i < numip; ++i ) {
159 ipindx(i) = pindx(i);
160 ipcoords(i,0) = pcoords(i,0);
161 ipcoords(i,1) = pcoords(i,1);
166 if(out.
get()) *out <<
"\nGenerating the element list ...\n";
169 global_node_id = localNodeOffset+1;
170 for(
int i_y = 0; i_y < local_ny; ++i_y ) {
171 for(
int i_x = 0; i_x < local_nx; ++i_x ) {
172 t(ele_k,0) = global_node_id;
173 t(ele_k,1) = global_node_id + (local_nx + 1);
174 t(ele_k,2) = global_node_id + 1;
175 if(out.
get() && dumpAll)
176 *out <<
"ele_k="<<ele_k<<
",("<<t(ele_k,0)<<
","<<t(ele_k,1)<<
","<<t(ele_k,2)<<
")\n";
178 t(ele_k,0) = global_node_id + 1;
179 t(ele_k,1) = global_node_id + (local_nx + 1);
180 t(ele_k,2) = global_node_id + (local_nx + 1) + 1;
181 if(out.
get() && dumpAll)
182 *out <<
"ele_k="<<ele_k<<
",("<<t(ele_k,0)<<
","<<t(ele_k,1)<<
","<<t(ele_k,2)<<
")\n";
196 if(out.
get()) *out <<
"\nGenerating the bottom edges ...\n";
198 global_node_id = localNodeOffset+1;
199 for(
int i_x = 0; i_x < local_nx; ++i_x ) {
200 e(edge_j,0) = global_node_id;
201 e(edge_j,1) = global_node_id + 1;
202 e(edge_j,2) = bndy_marker;
203 if(out.
get() && dumpAll)
204 *out <<
"edge_j="<<edge_j<<
",("<<e(edge_j,0)<<
","<<e(edge_j,1)<<
"),"<<e(edge_j,2)<<
"\n";
211 if(out.
get()) *out <<
"\nGenerating the left edges ...\n";
213 global_node_id = localNodeOffset+1;
214 for(
int i_y = 0; i_y < local_ny; ++i_y ) {
215 e(edge_j,0) = global_node_id;
216 e(edge_j,1) = global_node_id + (local_nx + 1);
217 e(edge_j,2) = bndy_marker;
218 if(out.
get() && dumpAll)
219 *out <<
"edge_j="<<edge_j<<
",("<<e(edge_j,0)<<
","<<e(edge_j,1)<<
"),"<<e(edge_j,2)<<
"\n";
221 global_node_id += (local_nx + 1);
225 if(out.
get()) *out <<
"\nGenerating the right edges ...\n";
227 global_node_id = localNodeOffset + 1 + local_nx;
228 for(
int i_y = 0; i_y < local_ny; ++i_y ) {
229 e(edge_j,0) = global_node_id;
230 e(edge_j,1) = global_node_id + (local_nx + 1);
231 e(edge_j,2) = bndy_marker;
232 if(out.
get() && dumpAll)
233 *out <<
"edge_j="<<edge_j<<
",("<<e(edge_j,0)<<
","<<e(edge_j,1)<<
"),"<<e(edge_j,2)<<
"\n";
235 global_node_id += (local_nx + 1);
238 if(procRank==numProc-1) {
240 if(out.
get()) *out <<
"\nGenerating the top edges ...\n";
242 global_node_id = localNodeOffset+1+(local_nx+1)*local_ny;
243 for(
int i_x = 0; i_x < local_nx; ++i_x ) {
244 e(edge_j,0) = global_node_id;
245 e(edge_j,1) = global_node_id + 1;
246 e(edge_j,2) = bndy_marker;
247 if(out.
get() && dumpAll)
248 *out <<
"edge_j="<<edge_j<<
",("<<e(edge_j,0)<<
","<<e(edge_j,1)<<
"),"<<e(edge_j,2)<<
"\n";
basic_OSTab< CharT, Traits > & incrTab(const int tabs=1)
void rect2DMeshGenerator(const int numProc, const int procRank, const double len_x, const double len_y, const int local_nx, const int local_ny, const int bndy_marker, Epetra_IntSerialDenseVector *ipindx_out, Epetra_SerialDenseMatrix *ipcoords_out, Epetra_IntSerialDenseVector *pindx_out, Epetra_SerialDenseMatrix *pcoords_out, Epetra_IntSerialDenseMatrix *t_out, Epetra_IntSerialDenseMatrix *e_out, std::ostream *out, const bool dumpAll)
Generate a simple rectangular 2D triangular mesh that is only partitioned between processors in the y...
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)