Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ExpList1DHomogeneous2D.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ExpList1DHomogeneous2D.cpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: An ExpList1D which is homogeneous in 2 directions and so
33 // uses much of the functionality from a ExpList2D and its daughters
34 //
35 ///////////////////////////////////////////////////////////////////////////////
36 
38 
39 namespace Nektar
40 {
41  namespace MultiRegions
42  {
43  // Forward declaration for typedefs
46  {
47  }
48 
49  // Constructor for ExpList1DHomogeneous2D to act as a Explist1D field
51  const LibUtilities::BasisKey &HomoBasis_y,
52  const LibUtilities::BasisKey &HomoBasis_z,
53  const NekDouble lhom_y,
54  const NekDouble lhom_z,
55  const bool useFFT,
56  const bool dealiasing,
57  const Array<OneD, ExpListSharedPtr> &points):
58  ExpListHomogeneous2D(pSession,HomoBasis_y,HomoBasis_z,lhom_y,lhom_z,useFFT,dealiasing)
59  {
60  int n,nel;
61 
62  int nyz = m_ny*m_nz;
63 
64  ASSERTL1(nyz == points.num_elements(),"Size of basis number of points and number of lines are not the same");
65 
66  for(n = 0; n < points.num_elements(); ++n)
67  {
68  m_lines[n] = points[n];
69  }
70 
71  // Setup Default optimisation information.
72  nel = 1;
75 
76  SetCoeffPhys();
77  }
78 
79  /**
80  * @param In ExpList1DHomogeneous2D object to copy.
81  */
84  {
85  for(int n = 0; n < m_lines.num_elements(); ++n)
86  {
87  m_lines[n] = In.m_lines[n];
88  }
89 
90  SetCoeffPhys();
91  }
92 
93  /**
94  * Destructor
95  */
97  {
98  }
99 
101  {
102  int i,n,cnt;
103  int ncoeffs_per_line = m_lines[0]->GetNcoeffs();
104  int npoints_per_line = m_lines[0]->GetTotPoints();
105 
106  int nyzlines = m_lines.num_elements();
107 
108  // Set total coefficients and points
109  m_ncoeffs = ncoeffs_per_line*nyzlines;
110  m_npoints = npoints_per_line*nyzlines;
111 
112  m_coeffs = Array<OneD, NekDouble> (m_ncoeffs);
113  m_phys = Array<OneD, NekDouble> (m_npoints);
114 
115  int nel = m_lines[0]->GetExpSize();
116  m_coeff_offset = Array<OneD,int>(nel*nyzlines);
117  m_phys_offset = Array<OneD,int>(nel*nyzlines);
118  m_offset_elmt_id = Array<OneD,int>(nel*nyzlines);
119  Array<OneD, NekDouble> tmparray;
120 
121  for(cnt = n = 0; n < nyzlines; ++n)
122  {
123  m_lines[n]->SetCoeffsArray(tmparray= m_coeffs + ncoeffs_per_line*n);
124  m_lines[n]->SetPhysArray(tmparray = m_phys + npoints_per_line*n);
125 
126  for(i = 0; i < nel; ++i)
127  {
128  m_coeff_offset[cnt] = m_lines[n]->GetCoeff_Offset(i) + n*ncoeffs_per_line;
129  m_phys_offset[cnt] = m_lines[n]->GetPhys_Offset(i) + n*npoints_per_line;
130  m_offset_elmt_id[cnt++] = m_lines[n]->GetOffset_Elmt_Id(i) + n*nel;
131  }
132  }
133  }
134 
136  Array<OneD, NekDouble> &xc0,
137  Array<OneD, NekDouble> &xc1,
138  Array<OneD, NekDouble> &xc2)
139  {
140  int n,m,j;
141  Array<OneD, NekDouble> tmp_xc;
142  int nylines = m_homogeneousBasis_y->GetNumPoints();
143  int nzlines = m_homogeneousBasis_z->GetNumPoints();
144  int npoints = 1;
145 
146  // Fill x-y-z-direction
147  Array<OneD, const NekDouble> pts_y = m_homogeneousBasis_y->GetZ();
148  Array<OneD, const NekDouble> pts_z = m_homogeneousBasis_z->GetZ();
149 
150  Array<OneD, NekDouble> x(npoints);
151  Array<OneD, NekDouble> y(nylines);
152  Array<OneD, NekDouble> z(nzlines);
153 
154  Vmath::Smul(nylines,m_lhom_y/2.0,pts_y,1,y,1);
155  Vmath::Sadd(nylines,m_lhom_y/2.0,y,1,y,1);
156 
157  Vmath::Smul(nzlines,m_lhom_z/2.0,pts_z,1,z,1);
158  Vmath::Sadd(nzlines,m_lhom_z/2.0,z,1,z,1);
159 
160  m_lines[0]->GetCoords(x);
161 
162 
163  for(m = 0; m < nzlines; ++m)
164  {
165  for(j = 0; j < nylines; ++j)
166  {
167  for(n = 0; n < npoints; ++n)
168  {
169  Vmath::Fill(1,x[n],tmp_xc = xc0 + n +(j*npoints) + (m*npoints*nylines), 1);
170  Vmath::Fill(1,y[j],tmp_xc = xc1 + n +(j*npoints) + (m*npoints*nylines), 1);
171  Vmath::Fill(1,z[m],tmp_xc = xc2 + n +(j*npoints) + (m*npoints*nylines), 1);
172  }
173  }
174  }
175  }
176 
177  /**
178  * The operation calls the 2D plane coordinates through the
179  * function ExpList#GetCoords and then evaluated the third
180  * coordinate using the member \a m_lhom
181  *
182  * @param coord_0 After calculation, the \f$x_1\f$ coordinate
183  * will be stored in this array.
184  *
185  * @param coord_1 After calculation, the \f$x_2\f$ coordinate
186  * will be stored in this array. This
187  * coordinate might be evaluated using the
188  * predefined value \a m_lhom
189  *
190  * @param coord_2 After calculation, the \f$x_3\f$ coordinate
191  * will be stored in this array. This
192  * coordinate is evaluated using the
193  * predefined value \a m_lhom
194  */
195  void ExpList1DHomogeneous2D::v_GetCoords(Array<OneD, NekDouble> &xc0,
196  Array<OneD, NekDouble> &xc1,
197  Array<OneD, NekDouble> &xc2)
198  {
199  int n,m,j;
200  Array<OneD, NekDouble> tmp_xc;
201  int npoints = 1;
202 
203  int nylines = m_homogeneousBasis_y->GetNumPoints();
204  int nzlines = m_homogeneousBasis_z->GetNumPoints();
205 
206  // Fill z-direction
207  Array<OneD, const NekDouble> pts_y = m_homogeneousBasis_y->GetZ();
208  Array<OneD, const NekDouble> pts_z = m_homogeneousBasis_z->GetZ();
209 
210  Array<OneD, NekDouble> x(npoints);
211  Array<OneD, NekDouble> y(nylines);
212  Array<OneD, NekDouble> z(nzlines);
213 
214  m_lines[0]->GetCoords(x);
215 
216  Vmath::Smul(nylines,m_lhom_y/2.0,pts_y,1,y,1);
217  Vmath::Sadd(nylines,m_lhom_y/2.0,y,1,y,1);
218 
219  Vmath::Smul(nzlines,m_lhom_z/2.0,pts_z,1,z,1);
220  Vmath::Sadd(nzlines,m_lhom_z/2.0,z,1,z,1);
221 
222  for(m = 0; m < nzlines; ++m)
223  {
224  for(j = 0; j < nylines; ++j)
225  {
226  for(n = 0; n < npoints; ++n)
227  {
228  Vmath::Fill(1,x[n],tmp_xc = xc0 + n +(j*npoints) + (m*npoints*nylines), 1);
229  Vmath::Fill(1,y[j],tmp_xc = xc1 + n +(j*npoints) + (m*npoints*nylines), 1);
230  Vmath::Fill(1,z[m],tmp_xc = xc2 + n +(j*npoints) + (m*npoints*nylines), 1);
231  }
232  }
233  }
234  }
235 
236 
237  /**
238  * Perform the 2D Forward transform of a set of points representing a plane of
239  * boundary conditions which are merely the collection of the boundary conditions
240  * coming from each 1D expansion.
241  * @param inarray The value of the BC on each point of the y-z homogeneous plane.
242  * @param outarray The value of the the coefficient of the 2D Fourier expansion
243  */
244  //void HomoFwdTrans2D(const Array<OneD, const NekDouble> &inarray, Array<OneD, NekDouble> &outarray)
245  //{
246 
247 
248  //}
249 
250  /**
251  * Write Tecplot Files Zone
252  * @param outfile Output file name.
253  * @param expansion Expansion that is considered
254  */
255  void ExpList1DHomogeneous2D::v_WriteTecplotZone(std::ofstream &outfile, int expansion)
256  {
257  int i,j;
258 
259  int nquad0 = 1;
260  int nquad1 = m_homogeneousBasis_y->GetNumPoints();
261  int nquad2 = m_homogeneousBasis_z->GetNumPoints();
262 
263  Array<OneD,NekDouble> coords[3];
264 
265  coords[0] = Array<OneD,NekDouble>(3*nquad0*nquad1*nquad2);
266  coords[1] = coords[0] + nquad0*nquad1*nquad2;
267  coords[2] = coords[1] + nquad0*nquad1*nquad2;
268 
269  GetCoords(expansion,coords[0],coords[1],coords[2]);
270 
271  outfile << "Zone, I=" << nquad1 << ", J=" << nquad0*nquad2
272  << ", F=Block" << std::endl;
273 
274  for(j = 0; j < nquad1; ++j)
275  {
276  for(i = 0; i < nquad2*GetCoordim(0)+1; ++i)
277  {
278  outfile << coords[j][i] << " ";
279  }
280  outfile << std::endl;
281  }
282  }
283 
284 
285  void ExpList1DHomogeneous2D::v_WriteVtkPieceHeader(std::ofstream &outfile, int expansion)
286  {
287  int i,j;
288 
289  int nquad0 = 1;
290  int nquad1 = m_homogeneousBasis_y->GetNumPoints();
291  int nquad2 = m_homogeneousBasis_z->GetNumPoints();
292 
293  int ntot = nquad0*nquad1*nquad2;
294  int ntotminus = (nquad0)*(nquad1-1)*(nquad2-1);
295 
296  Array<OneD,NekDouble> coords[3];
297  coords[0] = Array<OneD,NekDouble>(ntot);
298  coords[1] = Array<OneD,NekDouble>(ntot);
299  coords[2] = Array<OneD,NekDouble>(ntot);
300  GetCoords(expansion,coords[0],coords[1],coords[2]);
301 
302  outfile << " <Piece NumberOfPoints=\""
303  << ntot << "\" NumberOfCells=\""
304  << ntotminus << "\">" << endl;
305  outfile << " <Points>" << endl;
306  outfile << " <DataArray type=\"Float32\" "
307  << "NumberOfComponents=\"3\" format=\"ascii\">" << endl;
308  outfile << " ";
309  for (i = 0; i < ntot; ++i)
310  {
311  for (j = 0; j < 3; ++j)
312  {
313  outfile << coords[j][i] << " ";
314  }
315  outfile << endl;
316  }
317  outfile << endl;
318  outfile << " </DataArray>" << endl;
319  outfile << " </Points>" << endl;
320  outfile << " <Cells>" << endl;
321  outfile << " <DataArray type=\"Int32\" "
322  << "Name=\"connectivity\" format=\"ascii\">" << endl;
323  for (i = 0; i < nquad0; ++i)
324  {
325  for (j = 0; j < nquad1-1; ++j)
326  {
327  outfile << j*nquad0 + i << " "
328  << j*nquad0 + i + 1 << " "
329  << (j+1)*nquad0 + i + 1 << " "
330  << (j+1)*nquad0 + i << endl;
331  }
332  }
333  outfile << endl;
334  outfile << " </DataArray>" << endl;
335  outfile << " <DataArray type=\"Int32\" "
336  << "Name=\"offsets\" format=\"ascii\">" << endl;
337  for (i = 0; i < ntotminus; ++i)
338  {
339  outfile << i*4+4 << " ";
340  }
341  outfile << endl;
342  outfile << " </DataArray>" << endl;
343  outfile << " <DataArray type=\"UInt8\" "
344  << "Name=\"types\" format=\"ascii\">" << endl;
345  for (i = 0; i < ntotminus; ++i)
346  {
347  outfile << "9 ";
348  }
349  outfile << endl;
350  outfile << " </DataArray>" << endl;
351  outfile << " </Cells>" << endl;
352  outfile << " <PointData>" << endl;
353  }
354 
355 
356  } //end of namespace
357 } //end of namespace
358 
359 
360 /**
361 * $Log: v $
362 *
363 **/
364