Rocstar  1.0
Rocstar multiphysics simulation application
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Rocout_cgns.C
Go to the documentation of this file.
1 /* *******************************************************************
2  * Rocstar Simulation Suite *
3  * Copyright@2015, Illinois Rocstar LLC. All rights reserved. *
4  * *
5  * Illinois Rocstar LLC *
6  * Champaign, IL *
7  * www.illinoisrocstar.com *
8  * sales@illinoisrocstar.com *
9  * *
10  * License: See LICENSE file in top level of distribution package or *
11  * http://opensource.org/licenses/NCSA *
12  *********************************************************************/
13 /* *******************************************************************
14  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, *
15  * EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES *
16  * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND *
17  * NONINFRINGEMENT. IN NO EVENT SHALL THE CONTRIBUTORS OR *
18  * COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER *
19  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, *
20  * ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
21  * USE OR OTHER DEALINGS WITH THE SOFTWARE. *
22  *********************************************************************/
26 /* Author John Norris
27  * Initial date: Oct. 19, 2004
28  */
29 
30 #include <algorithm>
31 #include <cmath>
32 #include <cstring>
33 #include <fstream>
34 #include <iomanip>
35 #include <iostream>
36 #include <limits>
37 #include <sstream>
38 #include <string>
39 #include <vector>
40 
41 #include "cgnslib.h"
42 #include "Rocout.h"
43 #include "Element_accessors.h"
44 #include "Rocout_cgns.h"
45 #include <unistd.h>
46 
47 // MS
48 char cwd[FILENAME_MAX];
49 #define GetCurrentDir getcwd
50 // MS End
51 #ifndef DOXYGEN_SHOULD_SKIP_THIS
52 USE_COM_NAME_SPACE
53 #endif
54 
55 // #define DEBUG_DUMP_PREFIX "."
56 //#define DEBUG_MSG(x) std::cout << "ROCOUT DEBUG: " << __LINE__ << ": " << x << std::endl
57 #define DEBUG_MSG(x)
58 
68 #define CG_CHECK(routine, args) \
69 { \
70  int ier = routine args; \
71  if (ier != 0 && errorhandle != "ignore") { \
72  std::cerr << "Rocout::write_attribute: " #routine " (line " \
73  << __LINE__ << " in " << __FILE__ << ") failed: " \
74  << cg_get_error() << std::endl; \
75  if (errorhandle == "abort") { \
76  if (COMMPI_Initialized()) \
77  MPI_Abort(MPI_COMM_WORLD, 0); \
78  else \
79  abort(); \
80  } \
81  } \
82 }
83 
93 #define SwitchOnCOMDataType(dType, funcCall) \
94  switch (dType) { \
95  case COM_CHAR: \
96  case COM_CHARACTER: { \
97  typedef char COM_TT; \
98  funcCall; \
99  break; } \
100  case COM_INT: \
101  case COM_INTEGER: { \
102  typedef int COM_TT; \
103  funcCall; \
104  break; } \
105  case COM_FLOAT: \
106  case COM_REAL: { \
107  typedef float COM_TT; \
108  funcCall; \
109  break; } \
110  case COM_DOUBLE: \
111  case COM_DOUBLE_PRECISION: { \
112  typedef double COM_TT; \
113  funcCall; \
114  break; } \
115  }
116 
122 class AutoCloser
123 {
124  public:
125  inline AutoCloser(int fn, const std::string& eh)
126  : m_fn(fn), errorhandle(eh) {}
127  inline ~AutoCloser() { CG_CHECK(cg_close, (m_fn)); }
128 
129  private:
130  int m_fn;
131  const std::string& errorhandle;
132 };
133 
139 class AutoCDer
140 {
141  public:
142  inline AutoCDer()
143  : m_cwd(getcwd(NULL, 0)) {}
144  inline ~AutoCDer() { chdir(m_cwd); free(m_cwd); }
145 
146  private:
147  char* m_cwd;
148 };
149 
156 static DataType_t COM2CGNS(int dType)
157 {
158  switch (dType) {
159  case COM_CHAR:
160  case COM_CHARACTER:
161  return Character;
162 
163  case COM_INT:
164  case COM_INTEGER:
165  return Integer;
166 
167  case COM_FLOAT:
168  case COM_REAL:
169  return RealSingle;
170 
171  case COM_DOUBLE:
172  case COM_DOUBLE_PRECISION:
173  return RealDouble;
174 
175  default:
176  std::cerr << "ROCOUT ERROR: Cannot convert COM type " << dType
177  << " to a CGNS type." << std::endl;
178  }
179 
180  return DataTypeNull;
181 }
182 
201 static int WriteElements(int fn, int B, int Z, int eType, int nElems,
202  const int* elems, bool isStaggered, int length,
203  int offset, const std::string& label,
204  const std::string& errorhandle)
205 {
206  ElementType_t elemType;
207  int nNodes;
208 
209  // Determine the number of nodes per element and the CGNS element type.
210  switch (eType) {
211  case Connectivity::BAR2:
212  elemType = BAR_2;
213  nNodes = 2;
214  break;
215 
216  case Connectivity::BAR3:
217  elemType = BAR_3;
218  nNodes = 3;
219  break;
220 
221  case Connectivity::TRI3:
222  elemType = TRI_3;
223  nNodes = 3;
224  break;
225 
226  case Connectivity::TRI6:
227  elemType = TRI_6;
228  nNodes = 6;
229  break;
230 
231  case Connectivity::QUAD4:
232  elemType = QUAD_4;
233  nNodes = 4;
234  break;
235 
236  case Connectivity::QUAD8:
237  elemType = QUAD_8;
238  nNodes = 8;
239  break;
240 
241  case Connectivity::QUAD9:
242  elemType = QUAD_9;
243  nNodes = 9;
244  break;
245 
246  case Connectivity::TET4:
247  elemType = TETRA_4;
248  nNodes = 4;
249  break;
250 
251  case Connectivity::TET10:
252  elemType = TETRA_10;
253  nNodes = 10;
254  break;
255 
256  case Connectivity::PYRIMID5:
257  elemType = PYRA_5;
258  nNodes = 5;
259  break;
260 
261  case Connectivity::PYRIMID14:
262  elemType = PYRA_14;
263  nNodes = 14;
264  break;
265 
266  case Connectivity::PRISM6:
267  elemType = PENTA_6;
268  nNodes = 6;
269  break;
270 
271  case Connectivity::PRISM15:
272  elemType = PENTA_15;
273  nNodes = 15;
274  break;
275 
276  case Connectivity::PRISM18:
277  elemType = PENTA_18;
278  nNodes = 18;
279  break;
280 
281  case Connectivity::HEX8:
282  elemType = HEXA_8;
283  nNodes = 8;
284  break;
285 
286  case Connectivity::HEX20:
287  elemType = HEXA_20;
288  nNodes = 20;
289  break;
290 
291  case Connectivity::HEX27:
292  elemType = HEXA_27;
293  nNodes = 27;
294  break;
295  };
296 
297  // Allocate memory for the new connectivity table.
298  std::vector<int> buffer;
299  const int* pData;
300 
301  // Loop through the elements and get the indices of their nodes.
302  if (elems != NULL) {
303  buffer.resize(std::max(nElems * nNodes, 1));
304  int i, j, k = 0;
305  if (isStaggered) {
306  for (i=0; i<nElems; ++i)
307  for (j=0; j<nNodes; ++j,++k)
308  buffer[k] = elems[j*length+i];
309  pData = &buffer[0];
310  } else {
311  pData = elems;
312  }
313  } else {
314  buffer.resize(1);
315  pData = &buffer[0];
316  ++nElems;
317  }
318 
319  // Write out the Elements_t node.
320  int S = 0;
321  CG_CHECK(cg_section_write, (fn, B, Z, label.c_str(), elemType, offset,
322  offset + nElems - 1, 0, &(pData[0]), &S));
323 
324  return S;
325 }
326 
341 template <typename T>
342 void FindRange(int rank, const int* size, const int* rind, const T* pData,
343  std::string& range)
344 {
345  DEBUG_MSG("Entering FindRange<T>");
346  int i, j, k;
347  const T* t = pData;
348  T min = std::numeric_limits<T>::max();
349  T max = -min;
350  for (k=0; k<(rank>2?size[2]-rind[5]:1); ++k) {
351  for (j=0; j<(rank>1?size[1]:1); ++j) {
352  for (i=0; i<size[0]; ++i,++t) {
353  if (i >= std::min(rind[0], size[0] - 1)
354  && i < std::max(1, size[0] - rind[1])
355  && j >= std::min(rind[2], size[1] - 1)
356  && j < std::max(1, size[1] - rind[3])
357  && k >= std::min(rind[4], size[2] - 1)) {
358  if (min > *t)
359  min = *t;
360  if (max < *t)
361  max = *t;
362  }
363  }
364  }
365  }
366 
367  std::ostringstream sout;
368  sout << min << ", " << max;
369  range = sout.str();
370 }
371 
372 void FindRange(int rank, const int* size, const int* rind, const char* pData,
373  std::string& range)
374 {
375  DEBUG_MSG("Entering FindRange<char>");
376  int i, j, k;
377  const char* t = pData;
378  char min = std::numeric_limits<char>::max();
379  char max = -min;
380  for (k=0; k<(rank>2?size[2]-rind[5]:1); ++k) {
381  for (j=0; j<(rank>1?size[1]:1); ++j) {
382  for (i=0; i<size[0]; ++i,++t) {
383  if (i >= std::min(rind[0], size[0] - 1)
384  && i < std::max(1, size[0] - rind[1])
385  && j >= std::min(rind[2], size[1] - 1)
386  && j < std::max(1, size[1] - rind[3])
387  && k >= std::min(rind[4], size[2] - 1)) {
388  if (min > *t)
389  min = *t;
390  if (max < *t)
391  max = *t;
392  }
393  }
394  }
395  }
396 
397  std::ostringstream sout;
398  sout << (int)min << ", " << (int)max;
399  range = sout.str();
400 }
401 
417 template <typename T>
418 void FindMagnitudeRange(int physDim, int rank, const int* size,
419  const int* rind, const T** pData, std::string& range)
420 {
421  int i, j, k;
422  std::vector<T> buf;
423  const T* t[3] = { pData[0], pData[1], pData[2] };
424  if (physDim == 2) {
425  int s = size[0];
426  if (rank == 2)
427  s *= size[1];
428  buf.resize(s, (T)0);
429  pData[2] = &(buf[0]);
430  }
431 
432  T min = std::numeric_limits<T>::max();
433  T max = (T)0;
434  for (k=0; k<(rank>2?size[2]-rind[5]:1); ++k) {
435  for (j=0; j<(rank>1?size[1]:1); ++j) {
436  for (i=0; i<size[0]; ++i,++t[0],++t[1],++t[2]) {
437  if (i >= std::min(rind[0], size[0] - 1)
438  && i < std::max(1, size[0] - rind[1])
439  && j >= std::min(rind[2], size[1] - 1)
440  && j < std::max(1, size[1] - rind[3])
441  && k >= std::min(rind[4], size[2] - 1)) {
442  T val = *t[0] * *t[0] + *t[1] * *t[1] + *t[2] * *t[2];
443  if (val < min)
444  min = val;
445  if (val > max)
446  max = val;
447  }
448  }
449  }
450  }
451 
452  std::ostringstream sout;
453  sout << std::sqrt((double)min) << ", " << std::sqrt((double)max);
454  range = sout.str();
455 }
456 
472 template <typename T>
473 void FindTraceRange(int physDim, int rank, const int* size, const int* rind,
474  const T** pData, std::string& range)
475 {
476  int i, j, k;
477  std::vector<T> buf;
478  const T* t[3];
479  t[0] = pData[0];
480  if (physDim == 2) {
481  t[1] = pData[3];
482  int s = size[0];
483  if (rank == 2)
484  s *= size[1];
485  buf.resize(s, (T)0);
486  t[2] = &(buf[0]);
487  } else {
488  t[1] = pData[4];
489  t[2] = pData[8];
490  }
491 
492  T min = std::numeric_limits<T>::max();
493  T max = -min;
494  for (k=0; k<(rank>2?size[2]-rind[5]:1); ++k) {
495  for (j=0; j<(rank>1?size[1]:1); ++j) {
496  for (i=0; i<size[0]; ++i,++t[0],++t[1],++t[2]) {
497  if (i >= std::min(rind[0], size[0] - 1)
498  && i < std::max(1, size[0] - rind[1])
499  && j >= std::min(rind[2], size[1] - 1)
500  && j < std::max(1, size[1] - rind[3])
501  && k >= std::min(rind[4], size[2] - 1)) {
502  T val = *t[0] + *t[1] + *t[2];
503  if (val < min)
504  min = val;
505  if (val > max)
506  max = val;
507  }
508  }
509  }
510  }
511 
512  std::ostringstream sout;
513  sout << min << ", " << max;
514  range = sout.str();
515 }
516 
518 
527 float ExtractExponent(std::string& unit)
528 {
529  if (unit[0] != '^')
530  return 1.f;
531 
532  std::istringstream in(&unit[1]);
533  float e;
534  in >> e;
535 
536  if (in.eof())
537  unit.erase();
538  else
539  in >> unit;
540 
541  return e;
542 }
543 
544 void ModifyExponents(const std::string& unit, float* exponents, float dir)
545 {
546  std::string u = unit;
547  while (!u.empty()) {
548  float exp;
549  if (u.substr(0, 2) == "kg") { // Kilograms
550  u.erase(0, 2);
551  exp = ExtractExponent(u);
552  exponents[0] += exp * dir;
553  } else if (u[0] == 'm') { // Meters
554  u.erase(0, 1);
555  exp = ExtractExponent(u);
556  exponents[1] += exp * dir;
557  } else if (u[0] == 's') { // Seconds
558  u.erase(0, 1);
559  exp = ExtractExponent(u);
560  exponents[2] += exp * dir;
561  } else if (u[0] == 'K') { // Kelvin
562  u.erase(0, 1);
563  exp = ExtractExponent(u);
564  exponents[3] += exp * dir;
565  } else if (u[0] == 'J') { // Joules (kg m^2 / s^2)
566  u.erase(0, 1);
567  exp = ExtractExponent(u);
568  exponents[0] += exp * dir;
569  exponents[1] += 2 * exp * dir;
570  exponents[2] -= 2 * exp * dir;
571  } else if (u[0] == 'N') { // Newtons (kg m / s^2)
572  u.erase(0, 1);
573  exp = ExtractExponent(u);
574  exponents[0] += exp * dir;
575  exponents[1] += exp * dir;
576  exponents[2] -= 2 * exp * dir;
577  } else if (u.substr(0, 2) == "Pa") { // Pascals (kg / m s^2)
578  u.erase(0, 2);
579  exp = ExtractExponent(u);
580  exponents[0] += exp * dir;
581  exponents[1] -= exp * dir;
582  exponents[2] -= 2 * exp * dir;
583  } else {
584  std::cerr << "ROCOUT ERROR: Unable to parse unit '" << unit
585  << "': stuck at '" << u << "'." << std::endl;
586  break;
587  }
588  }
589 }
590 
591 template <typename T>
592 static int cg_array_core_write_internal(char const* name, DataType_t dType,
593  int rank, int* rind, int const* size,
594  T const* data)
595 {
596  int i, j, k, x = 0, total = 1;
597  int lMin[3] = { 0, 0, 0 }, lMax[3] = { 1, 1, 1 };
598  int lSize[3] = { 1, 1, 1 }, newSize[3] = { 1, 1, 1 };
599  for (i=0; i<rank; ++i) {
600  lMin[i] = rind[2*i];
601  lMax[i] = size[i] - rind[2*i+1];
602  lSize[i] = size[i];
603  newSize[i] = lMax[i] - lMin[i];
604  total *= newSize[i];
605  }
606  std::vector<T> buffer(total);
607 
608  // Extract subset.
609  for (k=lMin[2]; k<lMax[2]; ++k)
610  for (j=lMin[1]; j<lMax[1]; ++j)
611  for (i=lMin[0]; i<lMax[0]; ++i,++x)
612  buffer[x] = data[i+j*lSize[0]+k*lSize[0]*lSize[1]];
613 
614  return cg_array_write(name, dType, rank, newSize, &buffer[0]);
615 }
616 
617 static int cg_array_core_write(char const* name, DataType_t dType,
618  int rank, int* rind, int const* size,
619  void const* data)
620 {
621  // Check for a full write.
622  bool full = true;
623  int i;
624  for (i=0; i<rank && full; ++i)
625  full = (rind[2*i] == 0 && rind[2*i+1] == 0);
626  if (full)
627  return cg_array_write(name, dType, rank, size, data);
628 
629  switch (dType) {
630  case Character:
631  return cg_array_core_write_internal(name, dType, rank, rind, size,
632  static_cast<char const*>(data));
633  case Integer:
634  return cg_array_core_write_internal(name, dType, rank, rind, size,
635  static_cast<int const*>(data));
636  case RealSingle:
637  return cg_array_core_write_internal(name, dType, rank, rind, size,
638  static_cast<float const*>(data));
639  case RealDouble:
640  return cg_array_core_write_internal(name, dType, rank, rind, size,
641  static_cast<double const*>(data));
642  }
643 
644  return CG_ERROR;
645 }
646 
647 
648 static void cg_exponents_as_string_write(std::string unit,
649  const std::string& errorhandle)
650 {
651  // Eliminate spaces and parentheses.
652  std::string::size_type i = 0;
653  while (i < unit.length()) {
654  if (unit[i] == ' ' || unit[i] == '(' || unit[i] == ')')
655  unit.erase(i, 1);
656  else
657  ++i;
658  }
659 
660  // Normalize exponent operator.
661  i = unit.find("**");
662  while (i != std::string::npos) {
663  unit.replace(i, 2, "^");
664  i = unit.find("**", i);
665  }
666 
667  // Find the numerator and denominator.
668  std::string top, bottom;
669  i = unit.find('/');
670  if (i == std::string::npos) {
671  top = unit;
672  } else {
673  top = unit.substr(0, i);
674  bottom = unit.substr(i + 1);
675  }
676 
677  float exponents[5] = { 0.f, 0.f, 0.f, 0.f, 0.f };
678  if (!top.empty() && top != "1")
679  ModifyExponents(top, exponents, 1.f);
680  if (!bottom.empty() && bottom != "1")
681  ModifyExponents(bottom, exponents, -1.f);
682 
683  CG_CHECK(cg_exponents_write, (RealSingle, exponents));
684 }
686 
701 static int cg_base_find_or_create(int fn, const char* name, int cellDim,
702  int physDim, int* B,
703  const std::string& errorhandle)
704 {
705  char baseName[33];
706  int nBases, cDim, pDim;
707  CG_CHECK(cg_nbases, (fn, &nBases));
708 
709  // Check each existing Base_t node for a name match.
710  for (*B=1; *B<=nBases; ++(*B)) {
711  CG_CHECK(cg_base_read, (fn, *B, baseName, &cDim, &pDim));
712  if (std::strncmp(name, baseName, 32) == 0) {
713  // Confirm that the dimensions are correct.
714  if (cDim == cellDim && pDim == physDim)
715  return 0;
716  DEBUG_MSG("ERROR: Dimensions don't match (cellDim == " << cDim
717  << ", physDim == " << pDim << ")");
718  return 1;
719  }
720  }
721 
722  CG_CHECK(cg_base_write, (fn, name, cellDim, physDim, B));
723  return 0;
724 }
725 
741 static int cg_zone_find_or_create(int fn, int B, const char* name,
742  const int* sizes, ZoneType_t zType, int* Z,
743  const std::string& errorhandle)
744 {
745  char zoneName[33];
746  int nZones;
747  std::vector<int> sz1(&(sizes[0]), &(sizes[9]));
748  std::vector<int> sz2(9);
749  CG_CHECK(cg_nzones, (fn, B, &nZones));
750  for (*Z=1; *Z<=nZones; ++(*Z)) {
751  CG_CHECK(cg_zone_read, (fn, B, *Z, zoneName, &(sz2[0])));
752  if (std::strncmp(name, zoneName, 32) == 0) {
753  if (sz1 != sz2) {
754  // MS
755  //DEBUG_MSG("ERROR: Zone dimensions don't match");
756  //return 1;
757  // MS End
758  }
759 
760  ZoneType_t zt;
761  CG_CHECK(cg_zone_type, (fn, B, *Z, &zt));
762  if (zType == zt)
763  return 0;
764 
765  DEBUG_MSG("ERROR: Zone type doesn't match");
766  return 1;
767  }
768  }
769 
770  if (zType == Unstructured && sizes[2] > sizes[0])
771  std::cerr << "ROCOUT ERROR: Creating unstructured zone \"" << name
772  << "\" with invalid dimensions { " << sizes[0] << ", " << sizes[1]
773  << ", " << sizes[2] << " }" << std::endl;
774 
775  CG_CHECK(cg_zone_write, (fn, B, name, sizes, zType, Z));
776  return 0;
777 }
778 
789 static int cg_integral_find_or_create(const char* name, int* I,
790  const std::string& errorhandle)
791 {
792  char intName[33];
793  int nIntegrals;
794  CG_CHECK(cg_nintegrals, (&nIntegrals));
795  for (*I=1; *I<=nIntegrals; ++(*I)) {
796  CG_CHECK(cg_integral_read, (*I, intName));
797  if (std::strncmp(name, intName, 32) == 0)
798  return 0;
799  }
800 
801  CG_CHECK(cg_integral_write, (name));
802  return 0;
803 }
804 
819 static int cg_sol_find_or_create(int fn, int B, int Z, const char* name,
820  GridLocation_t location, int* S,
821  const std::string& errorhandle)
822 {
823  char solName[33];
824  GridLocation_t loc;
825  int nSols;
826  CG_CHECK(cg_nsols, (fn, B, Z, &nSols));
827  for (*S=1; *S<=nSols; ++(*S)) {
828  CG_CHECK(cg_sol_info, (fn, B, Z, *S, solName, &loc));
829  if (std::strncmp(name, solName, 32) == 0) {
830  return (loc == location ? 0 : 1);
831  }
832  }
833 
834  CG_CHECK(cg_sol_write, (fn, B, Z, name, location, S));
835  return 0;
836 }
837 
852 static int cg_array_info_by_name(const char* name, int* A, DataType_t* dType,
853  int* rank, int* size,
854  const std::string& errorhandle)
855 {
856  char arrayName[33];
857  int nArrays;
858  CG_CHECK(cg_narrays, (&nArrays));
859  //std::cout << __FILE__<< __LINE__ << std::endl;
860  //std::cout << " nArrays = " << nArrays << std::endl;
861  for (*A=1; *A<=nArrays; ++(*A)) {
862  CG_CHECK(cg_array_info, (*A, arrayName, dType, rank, size));
863  //std::cout << " arrayName = " << arrayName << std::endl;
864  if (std::strncmp(name, arrayName, 32) == 0)
865  return 0;
866  }
867  return 1;
868 }
869 
884 void write_attr_CGNS(const std::string& fname_in, const std::string& mfile,
885  const COM::Attribute* attr, const char* material,
886  const char* timelevel, int pane_id,
887  const std::string& ghosthandle,
888  const std::string& errorhandle, int mode)
889 {
890  /*
891  std::cout << " ------------------------------------------------------" << std::endl;
892  std::cout << " Starting to write \n Data File = " << fname_in << std::endl;
893  std::cout << " Mesh File = " << mfile << std::endl;
894  std::cout << " timelevel = " << timelevel << std::endl;
895  std::cout << " pane_id = " << pane_id << std::endl;
896  std::cout << " mode = " << mode << std::endl;
897  std::cout << " ------------------------------------------------------" << std::endl;
898  */
899 
900  std::string fname(fname_in);
901  bool writeGhost = (ghosthandle == "write");
902  const Window* w = attr->window();
903  COM_assertion(w != NULL);
904  const Pane& pane = w->pane(pane_id);
905 
906 
907  // Convert the timelevel from a string to a number.
908  double timeValue;
909  {
910  std::istringstream sin(timelevel);
911  sin >> timeValue;
912  }
913 
914  // Covert the timeValue back to a string. This "normalizes" the string.
915  std::string timeLevel;
916  {
917  std::ostringstream sout;
918  sout << timeValue;
919  timeLevel = sout.str();
920  }
921 
922  DEBUG_MSG("Writing to file '" << fname << "', mode == '"
923  << (mode ? "append'" : "write'"));
924  DEBUG_MSG("Using mesh file '" << mfile << "'");
925  AutoCDer autoCD;
926  std::string::size_type loc = fname.rfind('/');
927  // MS
928  GetCurrentDir(cwd, sizeof(cwd));
929  //std::cout << __FILE__ << __LINE__ << "\n" << cwd << std::endl;
930  std::string prefix = fname.substr(0, loc+1);
931  // original
932  /*
933  if (loc != std::string::npos) {
934  std::cout << __FILE__ << __LINE__ << " <<<<<<<<<<<<<<";
935  chdir(fname.substr(0, loc).c_str());
936  fname.erase(0, loc + 1);
937  }
938  */
939  // original
940  // MS
941 
942  // Open or create the file.
943  int fn;
944  if (mode == 0) {
945  CG_CHECK(cg_open, (fname.c_str(), MODE_WRITE, &fn));
946  CG_CHECK(cg_close, (fn));
947  }
948  CG_CHECK(cg_open, (fname.c_str(), MODE_MODIFY, &fn));
949 
950  // The file will be closed automagically when we exit this function.
951  AutoCloser autoCloser(fn, errorhandle);
952 
953  // Find or create the base (corresponds to window/material).
954  int i, B, nSteps = 0;
955  std::vector<double> times;
956  char buffer[33];
957  std::string label;
958  int cellDim = pane.dimension();
959  int physDim = pane.attribute(COM::COM_NC)->size_of_components();
960  //std::cout << __FILE__ << __LINE__ << std::endl;
961  //std::cout << " cell dim = " << cellDim << std::endl;
962  //std::cout << " phys dim = " << physDim << std::endl;
963  // MS
964  if (cellDim == 0)
965  cellDim = 3;
966  // MS End
967  CG_CHECK(cg_base_find_or_create, (fn, material, cellDim, physDim, &B,
968  errorhandle));
969  //std::cout << __FILE__ << __LINE__ << std::endl;
970 
971  // Write the default Roccom units at the top.
972  CG_CHECK(cg_goto, (fn, B, "end"));
973  CG_CHECK(cg_units_write, (Kilogram, Meter, Second, Kelvin, Degree));
974 
975  // Read the time values in the BaseIterativeData_t node.
976  // MS
977  if (timeValue == 0){
978  times.resize(1);
979  } else {
980  if (mode > 0 && cg_biter_read(fn, B, buffer, &nSteps) == 0) {
981  CG_CHECK(cg_goto, (fn, B, "BaseIterativeData_t", 1, "end"));
982  times.resize(nSteps + 1);
983  CG_CHECK(cg_array_read_as, (1, RealDouble, &(times[0])));
984  } else {
985  times.resize(1);
986  }
987  }
988  // MS End
989  /* Original
990  if (mode > 0 && cg_biter_read(fn, B, buffer, &nSteps) == 0) {
991  CG_CHECK(cg_goto, (fn, B, "BaseIterativeData_t", 1, "end"));
992  times.resize(nSteps + 1);
993  CG_CHECK(cg_array_read_as, (1, RealDouble, &(times[0])));
994  } else {
995  times.resize(1);
996  }
997  */
998 
999  // Search for the current time.
1000  int timeIndex;
1001  for (timeIndex=0; timeIndex<nSteps; ++timeIndex)
1002  if (times[timeIndex] == timeValue) {
1003  //std::cout << " timeIndex = " << timeIndex << std::endl;
1004  break;
1005  }
1006 
1007  // Update the time values in the BaseIterativeData_t node if necessary.
1008  if (timeIndex == nSteps) {
1009  //std::cout << __FILE__ << __LINE__ << std::endl;
1010  times[timeIndex] = timeValue;
1011  ++nSteps;
1012 
1013  // Write/rewrite the time values to the BaseIterativeData_t node.
1014  CG_CHECK(cg_biter_write, (fn, B, "TimeIterValues", nSteps));
1015  CG_CHECK(cg_goto, (fn, B, "BaseIterativeData_t", 1, "end"));
1016  CG_CHECK(cg_array_write, ("TimeValues", RealDouble, 1, &nSteps,
1017  &(times[0])));
1018  --nSteps;
1019  }
1020 
1021  // Find or create the zone (corresponds to pane/section). Use the pane
1022  // id as the zone name.
1023  std::string zoneName;
1024  {
1025  std::ostringstream sout;
1026  sout << std::setw(4) << std::setfill('0') << pane.id();
1027  zoneName = sout.str();
1028  }
1029 
1030  // Set up the sizes based on dimensionality and zone type.
1031  int Z, sizes[9] = { 0, 0, 0, 0, 0, 0, 0, 0, 0 };
1032  ZoneType_t zType;
1033  i = 0;
1034  if (pane.is_structured()) {
1035  //std::cout << __FILE__ << __LINE__ << std::endl;
1036  // Sizes should be core nodes/elements only.
1037  int ghost = pane.size_of_ghost_layers();
1038  sizes[i++] = pane.size_i() - 2 * ghost;
1039  sizes[i++] = pane.size_j() - 2 * ghost;
1040  if (physDim > 2 && cellDim > 2)
1041  sizes[i++] = pane.size_k() - 2 * ghost;
1042  sizes[i++] = sizes[0] - 1;
1043  sizes[i++] = sizes[1] - 1;
1044  if (physDim > 2 && cellDim > 2)
1045  sizes[i++] = sizes[2] - 1;
1046  zType = Structured;
1047  } else {
1048  //std::cout << __FILE__ << __LINE__ << std::endl;
1049  sizes[0] = pane.size_of_real_nodes();
1050  sizes[1] = pane.size_of_real_elements();
1051  zType = Unstructured;
1052  }
1053 
1054  //std::cout << __FILE__ << __LINE__ << std::endl;
1055  // MS
1056  if (sizes[0] == 0){
1057  int meshfn;
1058  char zName[33];
1059  //std::vector<int> sz2(9);
1060  if (!mfile.empty()) {
1061  CG_CHECK(cg_open,
1062  ((prefix + mfile).c_str(), MODE_READ, &meshfn));
1063  cg_zone_read(meshfn, 1, 1, zName, &(sizes[0]));
1064  CG_CHECK(cg_close, (meshfn));
1065  } else {
1066  sizes[0] = 1;
1067  }
1068  }
1069  // MS End
1070  CG_CHECK(cg_zone_find_or_create, (fn, B, zoneName.c_str(), sizes,
1071  zType, &Z, errorhandle));
1072  //std::cout << __FILE__ << __LINE__ << std::endl;
1073 
1074  // Create the name for the GridCoordinates_t node.
1075  std::string gridName("Grid");
1076  if (timeLevel.length() + gridName.length() > 32)
1077  gridName.erase(32 - timeLevel.length());
1078  gridName += timeLevel;
1079 
1080  // Create the name for the IntegralData_t node for window attributes.
1081  std::string winName("WinData");
1082  if (timeLevel.length() + winName.length() > 32)
1083  winName.erase(32 - timeLevel.length());
1084  winName += timeLevel;
1085 
1086  // Create the name for the IntegralData_t node for pane attributes.
1087  std::string paneName("PaneData");
1088  if (timeLevel.length() + paneName.length() > 32)
1089  paneName.erase(32 - timeLevel.length());
1090  paneName += timeLevel;
1091 
1092  // Create the name for the IntegralData_t node for conn attributes.
1093  std::string connName("ConnData");
1094  if (timeLevel.length() + connName.length() > 32)
1095  connName.erase(32 - timeLevel.length());
1096  connName += timeLevel;
1097 
1098  // Create the name for the FlowSolution_t node for element attributes.
1099  std::string elemName("ElemData");
1100  if (timeLevel.length() + elemName.length() > 32)
1101  elemName.erase(32 - timeLevel.length());
1102  elemName += timeLevel;
1103 
1104  // Create the name for the FlowSolution_t node for node attributes.
1105  std::string nodeName("NodeData");
1106  if (timeLevel.length() + nodeName.length() > 32)
1107  nodeName.erase(32 - timeLevel.length());
1108  nodeName += timeLevel;
1109 
1110  // Read and update the grid coordinate and solution names in the
1111  // ZoneIterativeData_t node.
1112  char (*gridNames)[32] = NULL;
1113  char (*nodeNames)[32] = NULL;
1114  int rank, size[3];
1115  if (mode > 0 && cg_ziter_read(fn, B, Z, buffer) == 0) {
1116  //std::cout << __FILE__ << __LINE__ << std::endl;
1117  //std::cout << " Going to " << "Base = " << B << " Zone = "<< Z << std::endl;
1118  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "ZoneIterativeData_t", 1, "end"));
1119  int numStr = nSteps;
1120  if (timeIndex >= nSteps)
1121  ++numStr;
1122 
1123  gridNames = (char(*)[32])new char[32*numStr];
1124  std::fill_n((char*)gridNames, 32 * numStr, '\0');
1125  nodeNames = (char(*)[32])new char[32*numStr];
1126  std::fill_n((char*)nodeNames, 32 * numStr, '\0');
1127 
1128  int A;
1129  DataType_t dataType;
1130  // Add the grid coordinates name to the GridCoordinatesPointers array,
1131  // unless it's already there.
1132  CG_CHECK(cg_array_info_by_name,
1133  ("GridCoordinatesPointers", &A, &dataType, &rank, size,
1134  errorhandle));
1135  COM_assertion_msg((dataType == Character && rank == 2
1136  && size[1] == numStr),
1137  "GridCoordinatesPointers in ZoneIterativeData is corrupt");
1138  CG_CHECK(cg_array_read, (A, gridNames));
1139  if (timeIndex < nSteps) {
1140  COM_assertion_msg(gridName == gridNames[timeIndex],
1141  "GridCoordinatesPointers in ZoneIterativeData is corrupt");
1142  } else {
1143  std::strncpy(gridNames[timeIndex], gridName.c_str(), 32);
1144  }
1145 
1146  // Add the flow solutions name to the FlowSolutionsPointers array,
1147  // unless it's already there.
1148  CG_CHECK(cg_array_info_by_name,
1149  ("FlowSolutionsPointers", &A, &dataType, &rank, size,
1150  errorhandle));
1151  COM_assertion_msg((dataType == Character && rank == 2
1152  && size[1] == numStr),
1153  "FlowSolutionsPointers in ZoneIterativeData is corrupt");
1154  CG_CHECK(cg_array_read, (A, nodeNames));
1155  if (timeIndex < nSteps) {
1156  COM_assertion_msg(nodeName == nodeNames[timeIndex],
1157  "FlowSolutionsPointers in ZoneIterativeData is corrupt");
1158  } else {
1159  std::strncpy(nodeNames[timeIndex], nodeName.c_str(), 32);
1160  }
1161  //std::cout << __FILE__ << __LINE__ << std::endl;
1162  } else {
1163  //std::cout << __FILE__ << __LINE__ << std::endl;
1164  //std::cout << "Writing ZoneIterativeData" << std::endl;
1165  cg_ziter_write(fn, B, Z, "ZoneIterativeData");
1166 
1167  gridNames = (char(*)[32])new char[32];
1168  std::fill_n((char*)gridNames, 32, '\0');
1169  std::strncpy(gridNames[timeIndex], gridName.c_str(), 32);
1170 
1171  nodeNames = (char(*)[32])new char[32];
1172  std::fill_n((char*)nodeNames, 32, '\0');
1173  std::strncpy(nodeNames[timeIndex], nodeName.c_str(), 32);
1174  }
1175  //std::cout << __FILE__ << __LINE__ << std::endl;
1176  //std::cout << " timeIndex = " << timeIndex
1177  // << " nSteps = " << nSteps << std::endl;
1178  if (timeIndex == nSteps) {
1179  //std::cout << " Writing under ZoneIterativeData "
1180  // << " timeIndex = " << timeIndex
1181  // << " nSteps = " << nSteps << std::endl;
1182  size[0] = 32;
1183  size[1] = ++nSteps;
1184 
1185  // Write/rewrite the grid coordinate and solution names in the
1186  // ZoneIterativeData_t node.
1187  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "ZoneIterativeData_t", 1, "end"));
1188  //std::cout << __FILE__ << __LINE__ << std::endl;
1189  CG_CHECK(cg_array_write,
1190  ("GridCoordinatesPointers", Character, 2, size, gridNames));
1191  //std::cout << __FILE__ << __LINE__ << std::endl;
1192  CG_CHECK(cg_array_write,
1193  ("FlowSolutionsPointers", Character, 2, size, nodeNames));
1194  }
1195  //std::cout << __FILE__ << __LINE__ << std::endl;
1196 
1197  delete[] (char*)gridNames;
1198  delete[] (char*)nodeNames;
1199 
1200  std::vector<char> buf;
1201  int mfn = fn, mB = B, mZ = Z, nGC = 0, G = 0, numItems, dType;
1202  int rind[6] = { 0, 0, 0, 0, 0, 0 };
1203  // Build a CGNS path to the current zone.
1204  std::string path("/"), range;
1205  path += material;
1206  path += '/' + zoneName + '/';
1207 
1208  // Write the grid coordinates and the connectivity tables, if any.
1209  /*
1210  std::cout << " COM::COM_MESH = " << COM::COM_MESH << "\n"
1211  << " COM::COM_PMESH = " << COM::COM_PMESH << "\n"
1212  << " COM::COM_CONN = " << COM::COM_CONN << "\n"
1213  << " COM::COM_NC1 = " << COM::COM_NC1 << "\n"
1214  << " COM::COM_NC3 = " << COM::COM_NC3 << "\n"
1215  << " COM::COM_NC = " << COM::COM_NC << "\n"
1216  << " COM::COM_ALL = " << COM::COM_ALL << "\n"
1217  << " attr->id()= " << attr->id() << "\n"
1218  << " attr->fullname()= " << attr->fullname() << "\n"
1219  << " attr->location()= " << attr->location() << "\n"
1220  << " attr->is_panel()= " << attr->is_panel() << "\n"
1221  << " attr->is_windowed()= " << attr->is_windowed() << "\n"
1222  << " attr->is_elemental()= " << attr->is_elemental() << "\n"
1223  << " attr->is_nodal()= " << attr->is_nodal() << "\n"
1224  << " attr->size_of_real_items()= " << attr->size_of_real_items() << "\n"
1225  << " attr->size_of_ghost_items()= " << attr->size_of_ghost_items() << std::endl;
1226  */
1227  if (attr->id() == COM::COM_MESH || attr->id() == COM::COM_PMESH
1228  || attr->id() == COM::COM_CONN
1229  || (attr->id() >= COM::COM_NC1 && attr->id() <= COM::COM_NC3)
1230  || attr->id() == COM::COM_NC || attr->id() == COM::COM_ALL) {
1231  //std::cout << __FILE__ << __LINE__ << std::endl;
1232  // Find or create the grid coordinates and element tables. If they're in a
1233  // different file, then we should create links from the main file to the
1234  // mesh file nodes.
1235  if (!mfile.empty()) {
1236  DEBUG_MSG("Writing to mesh file '" << mfile << "', mode == '"
1237  << (mode ? "append'" : "write'"));
1238  GetCurrentDir(cwd, sizeof(cwd));
1239  //std::cout << __FILE__ << __LINE__ << std::endl;
1240  CG_CHECK(cg_open,
1241  ((prefix + mfile).c_str(), mode ? MODE_MODIFY : MODE_WRITE, &mfn));
1242  //std::cout << __FILE__ << __LINE__ << std::endl;
1243 
1244  //std::cout << __FILE__ << __LINE__ << std::endl;
1245  CG_CHECK(cg_base_find_or_create, (mfn, material, cellDim, physDim, &mB,
1246  errorhandle));
1247  //std::cout << __FILE__ << __LINE__ << std::endl;
1248 
1249  CG_CHECK(cg_zone_find_or_create, (mfn, mB, zoneName.c_str(),
1250  sizes, zType, &mZ, errorhandle));
1251  }
1252  //std::cout << __FILE__ << __LINE__ << std::endl;
1253 
1254  // The GridCoordinates node is referenced by the GridCoordinatesPointers,
1255  // so we should create it even if we don't put any data in it.
1256  // The first grid written must always be names "GridCoordinates".
1257  // Since we are naming our data based on timeLevel, we'll need to use
1258  // a link for the first one.
1259  CG_CHECK(cg_ngrids, (mfn, mB, mZ, &nGC));
1260 
1261  if (nGC == 0) {
1262  CG_CHECK(cg_grid_write, (mfn, mB, mZ, "GridCoordinates", &G));
1263  if (mfile.empty()) {
1264  CG_CHECK(cg_goto, (mfn, mB, "Zone_t", mZ, "end"));
1265  CG_CHECK(cg_link_write, (gridName.c_str(), "",
1266  (path + "GridCoordinates").c_str()));
1267  } else {
1268  // Create links from the data file to the mesh file.
1269  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
1270  CG_CHECK(cg_ngrids, (fn, B, Z, &nGC));
1271  if (nGC == 0) {
1272  CG_CHECK(cg_link_write, ("GridCoordinates", (mfile).c_str(),
1273  (path + "GridCoordinates").c_str()));
1274  }
1275  CG_CHECK(cg_link_write, (gridName.c_str(), (mfile).c_str(),
1276  (path + "GridCoordinates").c_str()));
1277  }
1278  } else {
1279  if (mfile.empty()) {
1280  CG_CHECK(cg_grid_write, (mfn, mB, mZ, gridName.c_str(), &G));
1281  } else {
1282  // Create links from the data file to the mesh file.
1283  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
1284  CG_CHECK(cg_ngrids, (fn, B, Z, &nGC));
1285  if (nGC == 0) {
1286  CG_CHECK(cg_link_write, ("GridCoordinates", (mfile).c_str(),
1287  (path + "GridCoordinates").c_str()));
1288  }
1289  CG_CHECK(cg_link_write, (gridName.c_str(), (mfile).c_str(),
1290  (path + "GridCoordinates").c_str()));
1291  }
1292  }
1293  //std::cout << __FILE__ << __LINE__ << std::endl;
1294 
1295  if (attr->id() == COM::COM_MESH || attr->id() == COM::COM_PMESH
1296  || attr->id() == COM::COM_NC || attr->id() == COM::COM_ALL
1297  || (attr->id() >= COM::COM_NC1 && attr->id() <= COM::COM_NC3)) {
1298  //std::cout << __FILE__ << __LINE__ << std::endl;
1299  // Write the actual mesh data.
1300  if (G > 0) {
1301  CG_CHECK(cg_goto, (mfn, mB, "Zone_t", mZ, "GridCoordinates_t", G,
1302  "end"));
1303 
1304  if (pane.is_structured()) {
1305  rank = cellDim;
1306  size[0] = pane.size_i();
1307  size[1] = pane.size_j();
1308  size[2] = (physDim == 3 && cellDim == 3 ? pane.size_k(): 1);
1309  numItems = size[0] * size[1] * size[2];
1310  std::fill_n(rind, 2 * physDim, pane.size_of_ghost_layers());
1311  } else {
1312  rank = 1;
1313  size[0] = numItems = pane.size_of_nodes();
1314  std::fill_n(rind, 6, 0);
1315  rind[1] = pane.size_of_ghost_nodes();
1316  }
1317 
1318  // Write rind data.
1319  if (writeGhost)
1320  CG_CHECK(cg_rind_write, (rind));
1321 
1322  // std::string ranges[physDim];
1323  std::string ranges[3];
1324  int n, start = 0, finish = physDim;
1325  if (attr->id() >= COM::COM_NC1 && attr->id() <= COM::COM_NC3) {
1326  start = attr->id() - COM::COM_NC1;
1327  finish = start + 1;
1328  }
1329  for (n=start; n<finish; ++n) {
1330  const Attribute* nc = pane.attribute(COM::COM_NC1 + n);
1331  dType = nc->data_type();
1332 
1333  // Check for null/empty/strided data.
1334  const void* pData = nc->pointer();
1335  bool isNull = false;
1336  if (numItems <= 0) {
1337  buf.resize(COM::Attribute::get_sizeof(dType, 1), '\0');
1338  std::fill(buf.begin(), buf.end(), '\0');
1339  pData = &(buf[0]);
1340  size[0] = size[1] = size[2] = 1;
1341  ranges[n] = "EMPTY";
1342  } else if (pData == NULL) {
1343  buf.resize(COM::Attribute::get_sizeof(dType,
1344  std::max(numItems, 1)), '\0');
1345  std::fill(buf.begin(), buf.end(), '\0');
1346  pData = &(buf[0]);
1347  // size[0] = size[1] = size[2] = 1;
1348  isNull = true;
1349  ranges[n] = "NULL";
1350  } else {
1351  int stride = nc->stride();
1352  if (stride > 1) {
1353  int dSize = COM::Attribute::get_sizeof(dType, 1);
1354  buf.resize(dSize * std::max(numItems, 1));
1355  for (i=0; i<numItems; ++i)
1356  std::memcpy(&(buf[i*dSize]),
1357  &(((const char*)pData)[i*stride*dSize]), dSize);
1358  pData = &(buf[0]);
1359  }
1360 
1361 #ifdef DEBUG_DUMP_PREFIX
1362  {
1363  std::ofstream fout((DEBUG_DUMP_PREFIX + std::string(material) + ".nc" + '_' + timeLevel + ".cgns").c_str(), std::ios_base::app);
1364  switch (COM2CGNS(dType)) {
1365  case Character:
1366  for (i=0; i<numItems; ++i)
1367  fout << i << " : " << (int)((const char*)pData)[i]
1368  << '\n';
1369  break;
1370  case Integer:
1371  for (i=0; i<numItems; ++i)
1372  fout << i << " : " << ((const int*)pData)[i] << '\n';
1373  break;
1374  case RealSingle:
1375  for (i=0; i<numItems; ++i)
1376  fout << i << " : " << ((const float*)pData)[i] << '\n';
1377  break;
1378  case RealDouble:
1379  for (i=0; i<numItems; ++i)
1380  fout << i << " : " << ((const double*)pData)[i]
1381  << '\n';
1382  break;
1383  default:
1384  break;
1385  }
1386  fout << "###########################################\n";
1387  }
1388 #endif // DEBUG_DUMP_PREFIX
1389 
1390  SwitchOnCOMDataType(dType, FindRange(rank, size, rind,
1391  (const COM_TT*)pData,
1392  ranges[n]));
1393  }
1394  label = "Coordinate";
1395  label += (char)('X' + n);
1396  if (writeGhost) {
1397  CG_CHECK(cg_array_write, (label.c_str(), COM2CGNS(dType), rank,
1398  size, pData));
1399  } else {
1400  CG_CHECK(cg_array_core_write, (label.c_str(), COM2CGNS(dType), rank,
1401  rind, size, pData));
1402  }
1403  DEBUG_MSG("Writing attribute '" << (char)('x' + n) << "-nc', id == "
1404  << nc->id() << ", components == "
1405  << nc->size_of_components() << ", location == '"
1406  << nc->location() << '\'');
1407  }
1408 
1409  // Write dimensional exponents and ranges.
1410  for (n=start; n<finish; ++n) {
1411  CG_CHECK(cg_goto, (mfn, mB, "Zone_t", mZ, "GridCoordinates_t", G,
1412  "DataArray_t", n + 1, "end"));
1413  CG_CHECK(cg_descriptor_write, ("Range", ranges[n].c_str()));
1414  DEBUG_MSG("Writing descriptor 'Range': '" << ranges[n] << '\'');
1415  if (!attr->unit().empty()) {
1416  cg_exponents_as_string_write(attr->unit().c_str(), errorhandle);
1417  DEBUG_MSG("Writing descriptor 'Units': '" << attr->unit()
1418  << '\'');
1419  CG_CHECK(cg_descriptor_write, ("Units", attr->unit().c_str()));
1420  } else { // Assume meters.
1421  cg_exponents_as_string_write("m", errorhandle);
1422  DEBUG_MSG("Writing descriptor 'Units': 'm'");
1423  CG_CHECK(cg_descriptor_write, ("Units", "m"));
1424  }
1425  }
1426  }
1427  }
1428 
1429  if (pane.is_unstructured() &&
1430  (attr->id() == COM::COM_MESH || attr->id() == COM::COM_PMESH
1431  || attr->id() == COM::COM_CONN || attr->id() == COM::COM_ALL)) {
1432  // Now do the element tables.
1433  int nS;
1434  CG_CHECK(cg_nsections, (mfn, mB, mZ, &nS));
1435 
1436  if (nS == 0) {
1437  if (!mfile.empty())
1438  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
1439 
1440  std::vector<const Connectivity*> conns;
1441  pane.elements(conns);
1442 
1443  int S = 0;
1444  i = 0;
1445  std::vector<const Connectivity*>::const_iterator c;
1446  for (c=conns.begin(); c!=conns.end(); ++c) {
1447  int nItems = (*c)->size_of_items();
1448  int nGhost = (*c)->size_of_ghost_items();
1449 
1450  DEBUG_MSG("Writing attribute '" << (*c)->name() << "', nItems == "
1451  << nItems << ", elemType == " << (*c)->element_type());
1452  if (nItems == 0) {
1453  label = "Empty" + (*c)->name();
1454  S = WriteElements(mfn, mB, mZ, (*c)->element_type(), nItems,
1455  NULL, false, 0, (*c)->index_offset() + 1, label,
1456  errorhandle);
1457  } else if ((*c)->pointer() == NULL) {
1458  label = "Null" + (*c)->name();
1459  S = WriteElements(mfn, mB, mZ, (*c)->element_type(), nItems,
1460  NULL, false, 0, (*c)->index_offset() + 1, label,
1461  errorhandle);
1462  } else {
1463  // Write out a table of real nodes.
1464  label = (*c)->name();
1465  S = WriteElements(mfn, mB, mZ, (*c)->element_type(), nItems,
1466  (*c)->pointer(), ((*c)->stride() == 1),
1467  (*c)->capacity(), (*c)->index_offset() + 1,
1468  label, errorhandle);
1469  }
1470 
1471  // Write out a table of ghost nodes.
1472  if (nGhost > 0) {
1473  //std::cout << __FILE__ << __LINE__ << " Writing Gost node data " << std::endl;
1474  int elemRind[2] = { 0, nGhost };
1475  CG_CHECK(cg_goto, (mfn, mB, "Zone_t", mZ, "Elements_t", S, "end"));
1476  CG_CHECK(cg_rind_write, (elemRind));
1477  }
1478 
1479  // Link from the real file to the mesh file.
1480  if (!mfile.empty()) {
1481  /*
1482  std::cout << __FILE__ << __LINE__
1483  << " label =" << label.c_str()
1484  << std::endl;
1485  std::cout << __FILE__ << __LINE__
1486  << " mfile =" << mfile.c_str()
1487  << std::endl;
1488  std::cout << __FILE__ << __LINE__
1489  << " path + label =" << (path + label).c_str()
1490  << std::endl;
1491  */
1492  CG_CHECK(cg_link_write, (label.c_str(), (mfile).c_str(),
1493  (path + label).c_str()));
1494  }
1495 
1496  }
1497  }
1498  }
1499  } else if (!mfile.empty()) {
1500  // If we have a mesh file, we must link to the eternal GridCoordinates_t
1501  // and Elements_t nodes, even if attr->id() doesn't include mesh data.
1502  // First link the GridCoordinates_t node.
1503  //std::cout << __FILE__ << __LINE__ << std::endl;
1504  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
1505  CG_CHECK(cg_ngrids, (fn, B, Z, &nGC));
1506 
1507  if (nGC == 0) {
1508  CG_CHECK(cg_link_write, ("GridCoordinates", (mfile).c_str(),
1509  (path + "GridCoordinates").c_str()));
1510  }
1511  if (nGC == 0) {
1512  if (fname_in.compare(mfile)!=0){
1513  //std::cout << __FILE__ << __LINE__ << std::endl;
1514  CG_CHECK(cg_link_write, (gridName.c_str(), (mfile).c_str(),
1515  (path + "GridCoordinates").c_str()));
1516  } else {
1517  //std::cout << __FILE__ << __LINE__ << std::endl;
1518  CG_CHECK(cg_link_write, (gridName.c_str(), "",
1519  (path + "GridCoordinates").c_str()));
1520  }
1521 
1522  // Then link any Elements_t nodes.
1523  if (pane.is_unstructured()) {
1524  CG_CHECK(cg_open,
1525  ((prefix+mfile).c_str(), MODE_MODIFY, &mfn));
1526 
1527  CG_CHECK(cg_base_find_or_create, (mfn, material, cellDim, physDim, &mB,
1528  errorhandle));
1529 
1530  //std::cout << __FILE__ << __LINE__ << std::endl;
1531  CG_CHECK(cg_zone_find_or_create, (mfn, mB, zoneName.c_str(),
1532  sizes, zType, &mZ, errorhandle));
1533 
1534  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
1535 
1536  int S, nS;
1537  CG_CHECK(cg_nsections, (mfn, mB, mZ, &nS));
1538 
1539  char label[33];
1540  ElementType_t eType;
1541  int min, max, nBoundary, parent;
1542  for (S=1; S<=nS; ++S) {
1543  CG_CHECK(cg_section_read, (mfn, mB, mZ, S, label, &eType, &min, &max,
1544  &nBoundary, &parent));
1545 
1546  CG_CHECK(cg_link_write, (label, (mfile).c_str(), (path + label).c_str()));
1547  }
1548  }
1549 
1550  }
1551 
1552  }
1553  //std::cout << __FILE__ << __LINE__ << std::endl;
1554 
1555  if (mfn != fn){
1556  //std::cout << __FILE__ << __LINE__ << std::endl;
1557  CG_CHECK(cg_close, (mfn));
1558  }
1559 
1560  // This node is referenced in the FlowSolutionPointers, so we should make
1561  // sure it exists even if its left empty.
1562  int T;
1563  //std::cout << " nodeName = " << nodeName << std::endl;
1564  //std::cout << " Vertex = " << Vertex << std::endl;
1565  CG_CHECK(cg_sol_find_or_create, (fn, B, Z, nodeName.c_str(), Vertex, &T,
1566  errorhandle));
1567  //std::cout << __FILE__ << __LINE__ << std::endl;
1568 
1569  if (attr->id() == COM::COM_CONN || attr->id() == COM::COM_NC
1570  || (attr->id() >= COM::COM_NC1 && attr->id() <= COM::COM_NC3))
1571  return;
1572  //std::cout << __FILE__ << __LINE__ << std::endl;
1573 
1574  // Write out attributes. Window attributes should be stored in IntegralData
1575  // nodes under the base node. Pane, element, and node attributes should be
1576  // stored under the zone node, pane attrs in IntegralData nodes, element and
1577  // node attributes in FlowSolution nodes.
1578 
1579  std::vector<const Attribute*> attrs;
1580  if (attr->id() != COM::COM_MESH && attr->id() != COM::COM_PMESH
1581  && attr->id() != COM::COM_PCONN && attr->id() != COM::COM_RIDGES)
1582  pane.attributes(attrs);
1583  // MS Listing existing attributes in the pane
1584  /*
1585  std::cout << " //////////////////////////////////////////////////// " << std::endl;
1586  std::vector<const Attribute*> paneAttrs;
1587  pane.attributes(paneAttrs);
1588  std::vector<const Attribute*>::const_iterator ii;
1589  std::cout << " Listing existing attributes in the pane : " << std::endl;
1590  for (ii = paneAttrs.begin(); ii!=paneAttrs.end(); ii++)
1591  std::cout << " Attribute = " << (*ii)->fullname() << std::endl;
1592  std::cout << " //////////////////////////////////////////////////// " << std::endl;
1593  */
1594  // MS End
1595  //std::cout << __FILE__ << __LINE__ << std::endl;
1596 
1597  if (attr->id() == COM::COM_ALL || attr->id() == COM::COM_PMESH) {
1598  //std::cout << __FILE__ << __LINE__ << std::endl;
1599  attrs.insert(attrs.begin(), pane.attribute(COM::COM_RIDGES));
1600  attrs.insert(attrs.begin(), pane.attribute(COM::COM_PCONN));
1601  } else if (attr->id() == COM::COM_MESH) {
1602  //std::cout << __FILE__ << __LINE__ << std::endl;
1603  attrs.insert(attrs.begin(), pane.attribute(COM::COM_RIDGES));
1604  } else if (COM::COM_PCONN) {
1605  //std::cout << __FILE__ << __LINE__ << std::endl;
1606  attrs.insert(attrs.begin(), pane.attribute(COM::COM_PCONN));
1607  } else if (attr->id() == COM::COM_RIDGES) {
1608  //std::cout << __FILE__ << __LINE__ << std::endl;
1609  attrs.insert(attrs.begin(), pane.attribute(COM::COM_RIDGES));
1610  }
1611  //std::cout << "Number of attributes to write = " << attrs.size() << std::endl;
1612  //std::cout << __FILE__ << __LINE__ << std::endl;
1613 
1614  std::vector<const Attribute*>::const_iterator begin = attrs.begin();
1615  std::vector<const Attribute*>::const_iterator end = attrs.end();
1616  if ((attr->id() == COM::COM_MESH || attr->id() > COM::COM_ALL) && attr->id() != COM::COM_PMESH
1617  && attr->id() != COM::COM_ATTS && attr->id() != COM::COM_ALL) {
1618  while (begin != end && (*begin)->id() != attr->id())
1619  ++begin;
1620  COM_assertion_msg(begin != end,
1621  "Rocout::write_attribute: can't find attribute in pane");
1622  end = begin + 1;
1623  }
1624 
1625  std::vector<const Attribute*>::const_iterator a;
1626  // MS
1627  bool goToNextItem;
1628  char arrName[33];
1629  int tt2, tt1;
1630  DataType_t dT;
1631  // MS End
1632  for (a=begin; a!=end; ++a) {
1633  DEBUG_MSG("Writing attribute '" << (*a)->name() << "', id == "
1634  << (*a)->id() << ", components == " << (*a)->size_of_components()
1635  << ", location == '" << (*a)->location() << "\', datatype == "
1636  << (*a)->data_type());
1637  goToNextItem = false;
1638  int A, size[3], nComp = (*a)->size_of_components(), offset;
1639  const void* pData[9];
1640  dType = (*a)->data_type();
1641  switch ((*a)->location()) {
1642  case 'w':
1643  //std::cout << __FILE__ << __LINE__ << std::endl;
1644  size[0] = numItems = (*a)->size_of_items();
1645  size[1] = size[2] = 1;
1646  rind[0] = 0;
1647  rind[1] = (*a)->size_of_ghost_items();
1648 
1649  CG_CHECK(cg_goto, (fn, B, "end"));
1650  CG_CHECK(cg_integral_find_or_create, (winName.c_str(), &T,
1651  errorhandle));
1652 
1653  CG_CHECK(cg_goto, (fn, B, "IntegralData_t", T, "end"));
1654  CG_CHECK(cg_narrays, (&offset));
1655  //++offset;
1656  // MS
1657  // If the attribute we are trying to write is already written to
1658  // the file we just skip it.
1659  //std::cout << __FILE__ << __LINE__ << std::endl;
1660  //std::cout << " Zone = " << Z << " IntegralData_t = " << T << std::endl;
1661  for (int ii=1; ii<=offset; ii++)
1662  {
1663  cg_array_info(ii, arrName, &dT, &tt1, &tt2);
1664  //std::cout << "offset = " << ii <<" arrayName = " << arrName << std::endl;
1665  if ((std::string(arrName)).find(std::string((*a)->name())) != std::string::npos) {
1666  //std::cout << "The request for duplicate dataitem is ignored." << std::endl;
1667  goToNextItem = true;
1668  }
1669  }
1670  if (goToNextItem)
1671  break;
1672  ++offset;
1673 
1674  for (A=0; A<nComp; ++A) {
1675  CG_CHECK(cg_goto, (fn, B, "IntegralData_t", T, "end"));
1676 
1677  label = (*a)->name();
1678 
1679  const Attribute* pa;
1680  if (nComp == 1)
1681  pa = pane.attribute((*a)->id());
1682  else {
1683  pa = pane.attribute((*a)->id() + A + 1);
1684  std::ostringstream sout;
1685  sout << label << '#' << A + 1 << "of" << nComp;
1686  label = sout.str();
1687  }
1688 
1689 #ifdef DEBUG_DUMP_PREFIX
1690  std::ofstream fout((DEBUG_DUMP_PREFIX + std::string(material) + '.' + (*a)->name() + '_' + timeLevel + ".cgns").c_str(), std::ios_base::app);
1691 #endif // DEBUG_DUMP_PREFIX
1692  bool isNull = false;
1693  pData[A] = pa->pointer();
1694  if (numItems <= (writeGhost ? 0 : rind[1])) {
1695  buf.resize(COM::Attribute::get_sizeof(dType, 1), '\0');
1696  std::fill(buf.begin(), buf.end(), '\0');
1697  pData[A] = &(buf[0]);
1698  size[0] = 1;
1699  if (!writeGhost)
1700  rind[1] = 0;
1701  range = "EMPTY";
1702  } else if (pData[A] == NULL) {
1703  buf.resize(COM::Attribute::get_sizeof(dType,
1704  std::max(size[0], 1)), '\0');
1705  std::fill(buf.begin(), buf.end(), '\0');
1706  pData[A] = &(buf[0]);
1707  isNull = true;
1708  size[0] = 1;
1709  range = "NULL";
1710  } else {
1711  int stride = pa->stride();
1712  if (stride > 1) {
1713  int dSize = COM::Attribute::get_sizeof(dType, 1);
1714  buf.resize(dSize * std::max(numItems, 1));
1715  for (i=0; i<numItems; ++i)
1716  std::memcpy(&(buf[i*dSize]),
1717  &(((const char*)pData[A])[i*stride*dSize]), dSize);
1718  pData[A] = &(buf[0]);
1719  }
1720 
1721 #ifdef DEBUG_DUMP_PREFIX
1722  {
1723  switch (COM2CGNS(dType)) {
1724  case Character:
1725  for (i=0; i<numItems; ++i)
1726  fout << i << " : " << (int)((const char*)(pData[A]))[i]
1727  << '\n';
1728  break;
1729  case Integer:
1730  for (i=0; i<numItems; ++i)
1731  fout << i << " : " << ((const int*)(pData[A]))[i] << '\n';
1732  break;
1733  case RealSingle:
1734  for (i=0; i<numItems; ++i)
1735  fout << i << " : " << ((const float*)(pData[A]))[i] << '\n';
1736  break;
1737  case RealDouble:
1738  for (i=0; i<numItems; ++i)
1739  fout << i << " : " << ((const double*)(pData[A]))[i]
1740  << '\n';
1741  break;
1742  default:
1743  break;
1744  }
1745  fout << "###########################################\n";
1746  }
1747 #endif // DEBUG_DUMP_PREFIX
1748 
1749  // SwitchOnCOMDataType(pa->data_type(),
1750  SwitchOnCOMDataType(dType,
1751  FindRange(1, size, rind,
1752  (const COM_TT*)pData[A], range));
1753  }
1754  if (writeGhost) {
1755  CG_CHECK(cg_array_write, (label.c_str(), COM2CGNS(dType),
1756  1, size, pData[A]));
1757  } else {
1758  CG_CHECK(cg_array_core_write, (label.c_str(), COM2CGNS(dType),
1759  1, rind, size, pData[A]));
1760  }
1761 
1762  CG_CHECK(cg_goto, (fn, B, "IntegralData_t", T, "DataArray_t",
1763  A + offset, "end"));
1764 
1765  CG_CHECK(cg_descriptor_write, ("Range", range.c_str()));
1766  DEBUG_MSG("Writing descriptor 'Range': '" << range << '\'');
1767  if (!pa->unit().empty()) {
1768  cg_exponents_as_string_write(attr->unit().c_str(), errorhandle);
1769  CG_CHECK(cg_descriptor_write, ("Units", pa->unit().c_str()));
1770  DEBUG_MSG("Writing descriptor 'Units': '" << pa->unit() << '\'');
1771  }
1772 
1773  // We can't put a Rind_t node under a DataArray_t node.
1774  if (writeGhost) {
1775  std::ostringstream sout;
1776  sout << rind[1];
1777  CG_CHECK(cg_descriptor_write, ("Ghost", sout.str().c_str()));
1778  DEBUG_MSG("Writing descriptor 'Ghost': '" << sout.str() << '\'');
1779  }
1780  }
1781  break;
1782 
1783  case 'p':
1784  size[0] = numItems = (*a)->size_of_items();
1785  size[1] = size[2] = 1;
1786  rind[0] = 0;
1787  rind[1] = (*a)->size_of_ghost_items();
1788  //std::cout << __FILE__ << __LINE__ << std::endl;
1789  DEBUG_MSG("numItems == " << numItems << ", numGhostItems == "
1790  << rind[1]);
1791 
1792  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
1793  CG_CHECK(cg_integral_find_or_create, (paneName.c_str(), &T,
1794  errorhandle));
1795 
1796  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "IntegralData_t", T, "end"));
1797  CG_CHECK(cg_narrays, (&offset));
1798  //++offset;
1799  // MS
1800  // If the attribute we are trying to write is already written to
1801  // the file we just skip it.
1802  //std::cout << __FILE__ << __LINE__ << std::endl;
1803  //std::cout << " Zone = " << Z << " IntegralData_t = " << T << std::endl;
1804  for (int ii=1; ii<=offset; ii++)
1805  {
1806  cg_array_info(ii, arrName, &dT, &tt1, &tt2);
1807  //std::cout << "offset = " << ii <<" arrayName = " << arrName << std::endl;
1808  if ((std::string(arrName)).find(std::string((*a)->name())) != std::string::npos) {
1809  //std::cout << "The request for duplicate dataitem is ignored." << std::endl;
1810  goToNextItem = true;
1811  }
1812  }
1813  if (goToNextItem)
1814  break;
1815  ++offset;
1816  //if (offset == 4)
1817  // return;
1818  // MS End
1819 
1820  for (A=0; A<nComp; ++A) {
1821  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "IntegralData_t", T, "end"));
1822 
1823  label = (*a)->name();
1824 
1825  const Attribute* pa;
1826  if (nComp == 1)
1827  pa = pane.attribute((*a)->id());
1828  else {
1829  pa = pane.attribute((*a)->id() + A + 1);
1830  std::ostringstream sout;
1831  sout << label << '#' << A + 1 << "of" << nComp;
1832  label = sout.str();
1833  }
1834 
1835 #ifdef DEBUG_DUMP_PREFIX
1836  std::ofstream fout((DEBUG_DUMP_PREFIX + std::string(material) + '.' + (*a)->name() + '_' + timeLevel + ".cgns").c_str(), std::ios_base::app);
1837 #endif // DEBUG_DUMP_PREFIX
1838  bool isNull = false;
1839  pData[A] = pa->pointer();
1840  if (numItems <= (writeGhost ? 0 : rind[1])) {
1841  buf.resize(COM::Attribute::get_sizeof(dType, 1), '\0');
1842  std::fill(buf.begin(), buf.end(), '\0');
1843  pData[A] = &(buf[0]);
1844  size[0] = 1;
1845  if (!writeGhost)
1846  rind[1] = 0;
1847  range = "EMPTY";
1848  } else if (pData[A] == NULL) {
1849  buf.resize(COM::Attribute::get_sizeof(dType,
1850  std::max(size[0], 1)), '\0');
1851  std::fill(buf.begin(), buf.end(), '\0');
1852  pData[A] = &(buf[0]);
1853  isNull = true;
1854  size[0] = 1;
1855  range = "NULL";
1856  } else {
1857  int stride = pa->stride();
1858  if (stride > 1) {
1859  int dSize = COM::Attribute::get_sizeof(dType, 1);
1860  buf.resize(dSize * std::max(numItems, 1));
1861  for (i=0; i<numItems; ++i)
1862  std::memcpy(&(buf[i*dSize]),
1863  &(((const char*)pData[A])[i*stride*dSize]), dSize);
1864  pData[A] = &(buf[0]);
1865  }
1866 
1867 #ifdef DEBUG_DUMP_PREFIX
1868  {
1869  switch (COM2CGNS(dType)) {
1870  case Character:
1871  for (i=0; i<numItems; ++i)
1872  fout << i << " : " << (int)((const char*)(pData[A]))[i]
1873  << '\n';
1874  break;
1875  case Integer:
1876  for (i=0; i<numItems; ++i)
1877  fout << i << " : " << ((const int*)(pData[A]))[i] << '\n';
1878  break;
1879  case RealSingle:
1880  for (i=0; i<numItems; ++i)
1881  fout << i << " : " << ((const float*)(pData[A]))[i] << '\n';
1882  break;
1883  case RealDouble:
1884  for (i=0; i<numItems; ++i)
1885  fout << i << " : " << ((const double*)(pData[A]))[i]
1886  << '\n';
1887  break;
1888  default:
1889  break;
1890  }
1891  fout << "###########################################\n";
1892  }
1893 #endif // DEBUG_DUMP_PREFIX
1894 
1895  SwitchOnCOMDataType(dType,
1896  FindRange(1, size, rind,
1897  (const COM_TT*)pData[A], range));
1898  }
1899  if (writeGhost) {
1900  CG_CHECK(cg_array_write, (label.c_str(), COM2CGNS(dType),
1901  1, size, pData[A]));
1902  } else {
1903  CG_CHECK(cg_array_core_write, (label.c_str(), COM2CGNS(dType),
1904  1, rind, size, pData[A]));
1905  }
1906  //std::cout << __FILE__ << __LINE__
1907  // << " B = " << B << " Z = " << Z
1908  // << " T = " << T << " A = " << A
1909  // << " offset = " << offset
1910  // << std::endl;
1911  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "IntegralData_t", T,
1912  "DataArray_t", A + offset, "end"));
1913 
1914  CG_CHECK(cg_descriptor_write, ("Range", range.c_str()));
1915  DEBUG_MSG("Writing descriptor 'Range' == '" << range << '\'');
1916  if (!pa->unit().empty()) {
1917  cg_exponents_as_string_write(attr->unit().c_str(), errorhandle);
1918  DEBUG_MSG("Writing descriptor 'Units' == '" << pa->unit() << '\'');
1919  CG_CHECK(cg_descriptor_write, ("Units", pa->unit().c_str()));
1920  }
1921 
1922  // We can't put a Rind_t node under a DataArray_t node.
1923  if (writeGhost) {
1924  std::ostringstream sout;
1925  sout << rind[1];
1926  CG_CHECK(cg_descriptor_write, ("Ghost", sout.str().c_str()));
1927  DEBUG_MSG("Writing descriptor 'Ghost' == '" << sout.str() << '\'');
1928  }
1929  }
1930 
1931  if ((*a)->id() == COM::COM_PCONN) {
1932  // Create CGNS connectivity node.
1933  } else if ((*a)->id() == COM::COM_RIDGES) {
1934  DEBUG_MSG("Special COM_RIDGES handling: zType == "
1935  << (zType == Structured ? "Structured" : "Unstructured")
1936  << ", size_of_components() == "
1937  << (*a)->size_of_components() << ", size_of_items() == "
1938  << (*a)->size_of_items());
1939  if (zType == Unstructured && (*a)->size_of_components() == 2
1940  && (*a)->size_of_items() > 0) {
1941  // Create parallel base and zone with BAR_2 elements.
1942  std::string newName(material);
1943  newName += "_ridges";
1944  int RB;
1945  DEBUG_MSG("Creating base \"" << newName
1946  << "\", cellDim == 1, physDim == " << physDim);
1947  CG_CHECK(cg_base_find_or_create, (fn, newName.c_str(), 1, physDim,
1948  &RB, errorhandle));
1949 
1950  // Write the default Roccom units at the top.
1951  CG_CHECK(cg_goto, (fn, RB, "end"));
1952  CG_CHECK(cg_units_write, (Kilogram, Meter, Second, Kelvin, Degree));
1953 
1954  std::string rpath("/");
1955  rpath += material;
1956  rpath += '/';
1957  // Check to see if a BaseIterativeData node exists.
1958  if (cg_goto(fn, RB, "BaseIterativeData_t", 1, "end") != 0) {
1959  // Link to the main BaseIterativeData_t node.
1960  DEBUG_MSG("Linking TimeIterValues to "
1961  << rpath + "TimeIterValues");
1962  CG_CHECK(cg_goto, (fn, RB, "end"));
1963  CG_CHECK(cg_link_write, ("TimeIterValues", "",
1964  (rpath + "TimeIterValues").c_str()));
1965  }
1966 
1967  // Find or create the zone (corresponds to pane/section). Use the
1968  // pane id as the zone name. Set up the sizes based on
1969  // dimensionality and zone type.
1970  int RZ, rsizes[9] = { 0, 0, 0, 0, 0, 0, 0, 0, 0 };
1971  rsizes[0] = sizes[0];
1972  rsizes[1] = (*a)->size_of_real_items();
1973  DEBUG_MSG("Creating zone \"" << zoneName << "\", size == [ "
1974  << rsizes[0] << ", " << rsizes[1] << " ]");
1975  //std::cout << __FILE__ << " line : " << __LINE__ << std::endl;
1976  CG_CHECK(cg_zone_find_or_create, (fn, RB, zoneName.c_str(), rsizes,
1977  Unstructured, &RZ, errorhandle));
1978 
1979  // Check to see if a ZoneIterativeData node exists.
1980  if (cg_goto(fn, RB, "Zone_t", RZ,
1981  "ZoneIterativeData_t", 1, "end") != 0) {
1982  CG_CHECK(cg_goto, (fn, RB, "Zone_t", RZ, "end"));
1983  rpath += zoneName;
1984  rpath += '/';
1985 
1986  // Link to the main ZoneIterativeData_t node.
1987  DEBUG_MSG("Linking ZoneIterativeData to "
1988  << rpath + "ZoneIterativeData");
1989  CG_CHECK(cg_link_write, ("ZoneIterativeData", "",
1990  (rpath + "ZoneIterativeData").c_str()));
1991 
1992  // Link to the main GridCoordinates_t node.
1993  DEBUG_MSG("Linking GridCoordinates to "
1994  << rpath + "GridCoordinates");
1995  CG_CHECK(cg_link_write, ("GridCoordinates", "",
1996  (rpath + "GridCoordinates").c_str()));
1997  DEBUG_MSG("Linking " << gridName << " to " << rpath + gridName);
1998  CG_CHECK(cg_link_write, (gridName.c_str(), "",
1999  (rpath + gridName).c_str()));
2000 
2001  // Link to the FlowSolution_t node with nodal data.
2002  DEBUG_MSG("Linking " << nodeName << " to "
2003  << rpath + nodeName);
2004  CG_CHECK(cg_link_write, (nodeName.c_str(), "",
2005  (rpath + nodeName).c_str()));
2006  }
2007 
2008  const Attribute* pa[2];
2009  pa[0] = pane.attribute((*a)->id() + 1);
2010  pa[1] = pane.attribute((*a)->id() + 2);
2011  DEBUG_MSG("Retreived the 'ridges' components: '" << pa[0]->name()
2012  << "' and '" << pa[1]->name() << '\'');
2013 
2014  // Now do the BAR_2 element tables.
2015  int RS = 0;
2016  int nItems = pa[0]->size_of_items();
2017  int nGhost = pa[0]->size_of_ghost_items();
2018  DEBUG_MSG("Writing " << nItems << " BAR_2 elements (with "
2019  << nGhost << " ghost elements)");
2020  COM_assertion(nItems == pa[1]->size_of_items());
2021  COM_assertion(nGhost == pa[1]->size_of_ghost_items());
2022  COM_assertion((pa[0]->pointer()==NULL) == (pa[1]->pointer()==NULL));
2023 
2024  if (!writeGhost)
2025  nItems = pa[0]->size_of_real_items();
2026 
2027  // Write out a table of real nodes.
2028  label = (*a)->name();
2029  std::vector<int> bar(2 * nItems);
2030  for (i=0; i<2; ++i) {
2031  int* ptr = (int*)pa[i]->pointer();
2032  int j, stride = std::max(pa[i]->stride(), 1);
2033  for (j=0; j<nItems; ++j)
2034  bar[2*j+i] = ptr[j*stride];
2035  }
2036  RS = WriteElements(mfn, RB, RZ, Connectivity::BAR2, nItems,
2037  &bar[0], false, 0, 1, label, errorhandle);
2038 
2039  // Write out a table of ghost elements.
2040  if (writeGhost) {
2041  int elemRind[2] = { 0, nGhost };
2042  CG_CHECK(cg_goto,
2043  (fn, RB, "Zone_t", RZ, "Elements_t", RS, "end"));
2044  CG_CHECK(cg_rind_write, (elemRind));
2045  }
2046  }
2047  } else if ((*a)->name() == "bcflag") {
2048  // Create CGNS boundary conditions node.
2049  // Boundary conditions might have other names, too.
2050  }
2051  break;
2052 
2053  case 'c':
2054  size[0] = numItems = (*a)->size_of_items();
2055  size[1] = size[2] = 1;
2056  rind[0] = 0;
2057  rind[1] = (*a)->size_of_ghost_items();
2058 
2059  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "end"));
2060  CG_CHECK(cg_integral_find_or_create, (connName.c_str(), &T,
2061  errorhandle));
2062 
2063  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "IntegralData_t", T, "end"));
2064  CG_CHECK(cg_narrays, (&offset));
2065  ++offset;
2066 
2067  for (A=0; A<nComp; ++A) {
2068  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "IntegralData_t", T, "end"));
2069 
2070  label = (*a)->name();
2071 
2072  const Attribute* pa;
2073  if (nComp == 1)
2074  pa = pane.attribute((*a)->id());
2075  else {
2076  pa = pane.attribute((*a)->id() + A + 1);
2077  std::ostringstream sout;
2078  sout << label << '#' << A + 1 << "of" << nComp;
2079  label = sout.str();
2080  }
2081 
2082 #ifdef DEBUG_DUMP_PREFIX
2083  std::ofstream fout((DEBUG_DUMP_PREFIX + std::string(material) + '.' + (*a)->name() + '_' + timeLevel + ".cgns").c_str(), std::ios_base::app);
2084 #endif // DEBUG_DUMP_PREFIX
2085  bool isNull = false;
2086  pData[A] = pa->pointer();
2087  if (numItems <= (writeGhost ? 0 : rind[1])) {
2088  buf.resize(COM::Attribute::get_sizeof(dType, 1), '\0');
2089  std::fill(buf.begin(), buf.end(), '\0');
2090  pData[A] = &(buf[0]);
2091  size[0] = 1;
2092  if (!writeGhost)
2093  rind[1] = 0;
2094  range = "EMPTY";
2095  } else if (pData[A] == NULL) {
2096  buf.resize(COM::Attribute::get_sizeof(dType,
2097  std::max(size[0], 1)), '\0');
2098  std::fill(buf.begin(), buf.end(), '\0');
2099  pData[A] = &(buf[0]);
2100  isNull = true;
2101  size[0] = 1;
2102  range = "NULL";
2103  } else {
2104  int stride = pa->stride();
2105  if (stride > 1) {
2106  int dSize = COM::Attribute::get_sizeof(dType, 1);
2107  buf.resize(dSize * std::max(numItems, 1));
2108  for (i=0; i<numItems; ++i)
2109  std::memcpy(&(buf[i*dSize]),
2110  &(((const char*)pData[A])[i*stride*dSize]), dSize);
2111  pData[A] = &(buf[0]);
2112  }
2113 
2114 #ifdef DEBUG_DUMP_PREFIX
2115  {
2116  switch (COM2CGNS(dType)) {
2117  case Character:
2118  for (i=0; i<numItems; ++i)
2119  fout << i << " : " << (int)((const char*)(pData[A]))[i]
2120  << '\n';
2121  break;
2122  case Integer:
2123  for (i=0; i<numItems; ++i)
2124  fout << i << " : " << ((const int*)(pData[A]))[i] << '\n';
2125  break;
2126  case RealSingle:
2127  for (i=0; i<numItems; ++i)
2128  fout << i << " : " << ((const float*)(pData[A]))[i] << '\n';
2129  break;
2130  case RealDouble:
2131  for (i=0; i<numItems; ++i)
2132  fout << i << " : " << ((const double*)(pData[A]))[i]
2133  << '\n';
2134  break;
2135  default:
2136  break;
2137  }
2138  fout << "###########################################\n";
2139  }
2140 #endif // DEBUG_DUMP_PREFIX
2141 
2142  SwitchOnCOMDataType(dType,
2143  FindRange(1, size, rind,
2144  (const COM_TT*)pData[A], range));
2145  }
2146  if (writeGhost) {
2147  CG_CHECK(cg_array_write, (label.c_str(), COM2CGNS(dType),
2148  1, size, pData[A]));
2149  } else {
2150  CG_CHECK(cg_array_core_write, (label.c_str(), COM2CGNS(dType),
2151  1, rind, size, pData[A]));
2152  }
2153 
2154  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "IntegralData_t", T,
2155  "DataArray_t", A + offset, "end"));
2156 
2157  CG_CHECK(cg_descriptor_write, ("Range", range.c_str()));
2158  DEBUG_MSG("Writing descriptor 'Range': '" << range << '\'');
2159  if (!pa->unit().empty()) {
2160  cg_exponents_as_string_write(attr->unit().c_str(), errorhandle);
2161  CG_CHECK(cg_descriptor_write, ("Units", pa->unit().c_str()));
2162  DEBUG_MSG("Writing descriptor 'Units': '" << pa->unit() << '\'');
2163  }
2164 
2165  // We can't put a Rind_t node under a DataArray_t node.
2166  if (writeGhost) {
2167  std::ostringstream sout;
2168  sout << rind[1];
2169  CG_CHECK(cg_descriptor_write, ("Ghost", sout.str().c_str()));
2170  DEBUG_MSG("Writing descriptor 'Ghost': '" << sout.str() << '\'');
2171  }
2172  }
2173  break;
2174 
2175  case 'n':
2176  if (pane.is_structured()) {
2177  rank = cellDim;
2178  size[0] = pane.size_i();
2179  if (cellDim > 1) {
2180  size[1] = pane.size_j();
2181  size[2] = cellDim > 2 ? pane.size_k() : 1;
2182  } else {
2183  size[1] = size[2] = 1;
2184  }
2185  numItems = size[0] * size[1] * size[2];
2186  std::fill_n(rind, 2 * physDim, pane.size_of_ghost_layers());
2187  } else {
2188  rank = 1;
2189  size[0] = numItems = (*a)->size_of_items();
2190  size[1] = size[2] = 1;
2191  std::fill_n(rind, 6, 0);
2192  rind[1] = (*a)->size_of_ghost_items();
2193  }
2194 
2195  CG_CHECK(cg_sol_find_or_create, (fn, B, Z, nodeName.c_str(), Vertex, &T,
2196  errorhandle));
2197 
2198  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T, "end"));
2199  CG_CHECK(cg_gridlocation_write, (Vertex));
2200 
2201  // Write rind data.
2202  if (writeGhost) {
2203  CG_CHECK(cg_rind_write, (rind));
2204  DEBUG_MSG("Ghost == " << rind[1]);
2205  }
2206 
2207  CG_CHECK(cg_narrays, (&offset));
2208  ++offset;
2209 
2210  for (A=0; A<nComp; ++A) {
2211  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T, "end"));
2212 
2213  label = (*a)->name();
2214 
2215  const Attribute* pa;
2216  if (nComp == 1)
2217  pa = pane.attribute((*a)->id());
2218  else {
2219  pa = pane.attribute((*a)->id() + A + 1);
2220  if (nComp == physDim) {
2221  label += (char)('X' + (char)A);
2222  } else if (nComp == physDim * physDim) {
2223  label += (char)('X' + (char)(A / physDim));
2224  label += (char)('X' + (char)(A % physDim));
2225  } else {
2226  std::ostringstream sout;
2227  sout << label << '#' << A + 1 << "of" << nComp;
2228  label = sout.str();
2229  }
2230  }
2231 
2232 #ifdef DEBUG_DUMP_PREFIX
2233  std::ofstream fout((DEBUG_DUMP_PREFIX + std::string(material) + '.' + (*a)->name() + '_' + timeLevel + ".cgns").c_str(), std::ios_base::app);
2234 #endif // DEBUG_DUMP_PREFIX
2235  bool isNull = false;
2236  pData[A] = pa->pointer();
2237  if (numItems <= 0) {
2238  buf.resize(COM::Attribute::get_sizeof(dType, 1), '\0');
2239  std::fill(buf.begin(), buf.end(), '\0');
2240  pData[A] = &(buf[0]);
2241  size[0] = size[1] = size[2] = 1;
2242  range = "EMPTY";
2243  } else if (pData[A] == NULL) {
2244  buf.resize(COM::Attribute::get_sizeof(dType, std::max(numItems, 1)),
2245  '\0');
2246  std::fill(buf.begin(), buf.end(), '\0');
2247  pData[A] = &(buf[0]);
2248  isNull = true;
2249  // size[0] = size[1] = size[2] = 1;
2250  range = "NULL";
2251  } else {
2252  int stride = pa->stride();
2253  if (stride > 1) {
2254  int dSize = COM::Attribute::get_sizeof(dType, 1);
2255  buf.resize(dSize * std::max(numItems, 1));
2256  for (i=0; i<numItems; ++i) {
2257  std::memcpy(&(buf[i*dSize]),
2258  &(((const char*)pData[A])[i*stride*dSize]), dSize);
2259  }
2260  pData[A] = &(buf[0]);
2261  }
2262 #ifdef DEBUG_DUMP_PREFIX
2263  {
2264  switch (COM2CGNS(dType)) {
2265  case Character:
2266  for (i=0; i<numItems; ++i)
2267  fout << i << " : " << (int)((const char*)(pData[A]))[i]
2268  << '\n';
2269  break;
2270  case Integer:
2271  for (i=0; i<numItems; ++i)
2272  fout << i << " : " << ((const int*)(pData[A]))[i] << '\n';
2273  break;
2274  case RealSingle:
2275  for (i=0; i<numItems; ++i)
2276  fout << i << " : " << ((const float*)(pData[A]))[i] << '\n';
2277  break;
2278  case RealDouble:
2279  for (i=0; i<numItems; ++i)
2280  fout << i << " : " << ((const double*)(pData[A]))[i]
2281  << '\n';
2282  break;
2283  default:
2284  break;
2285  }
2286  fout << "###########################################\n";
2287  }
2288 #endif // DEBUG_DUMP_PREFIX
2289 
2290  SwitchOnCOMDataType(dType,
2291  FindRange(rank, size, rind,
2292  (const COM_TT*)pData[A], range));
2293  }
2294  if (writeGhost) {
2295  DEBUG_MSG("Calling cg_array_write( name == '" << label
2296  << "', dataType == "
2297  << (COM2CGNS(dType) == RealSingle ? "float"
2298  : (COM2CGNS(dType) == RealDouble ? "double"
2299  : (COM2CGNS(dType) == Character ? "char" : "int?")))
2300  << ", rank == " << rank << ", size[] = { " << size[0]
2301  << ", " << size[1] << ", " << size[2] << " } )");
2302  CG_CHECK(cg_array_write, (label.c_str(), COM2CGNS(dType),
2303  rank, size, pData[A]));
2304  } else {
2305  DEBUG_MSG("Calling cg_array_core_write( name == '" << label
2306  << "', dataType == "
2307  << (COM2CGNS(dType) == RealSingle ? "float"
2308  : (COM2CGNS(dType) == RealDouble ? "double"
2309  : (COM2CGNS(dType) == Character ? "char" : "int?")))
2310  << ", rank == " << rank << ", rind[] = { " << rind[0]
2311  << ", " << rind[1] << ", " << rind[2] << ", " << rind[3]
2312  << ", " << rind[4] << ", " << rind[5] << " }, size[] = { "
2313  << size[0] << ", " << size[1] << ", " << size[2]
2314  << " } )");
2315  CG_CHECK(cg_array_core_write, (label.c_str(), COM2CGNS(dType),
2316  rank, rind, size, pData[A]));
2317  }
2318 
2319  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T,
2320  "DataArray_t", A + offset, "end"));
2321 
2322  DEBUG_MSG("Calling cg_descriptor_write( name == 'Range', "
2323  << "value == '" << range << "' )");
2324  CG_CHECK(cg_descriptor_write, ("Range", range.c_str()));
2325  if (!pa->unit().empty()) {
2326  cg_exponents_as_string_write(attr->unit().c_str(), errorhandle);
2327  DEBUG_MSG("Calling cg_descriptor_write( name == 'Units', "
2328  << "value == '" << pa->unit() << "' )");
2329  CG_CHECK(cg_descriptor_write, ("Units", pa->unit().c_str()));
2330  }
2331  }
2332 
2333  // Vectors and tensors need a MagnitudeRange or TraceRange
2334  // descriptor under the first DataArray_t node.
2335  if (nComp == 3) {
2336  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T,
2337  "DataArray_t", offset, "end"));
2338  SwitchOnCOMDataType(dType,
2339  FindMagnitudeRange(physDim, rank, size, rind,
2340  (const COM_TT**)pData,
2341  range));
2342  CG_CHECK(cg_descriptor_write, ("MagnitudeRange", range.c_str()));
2343  } else if (nComp == 9) {
2344  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T,
2345  "DataArray_t", offset, "end"));
2346  SwitchOnCOMDataType(dType,
2347  FindTraceRange(physDim, rank, size, rind,
2348  (const COM_TT**)pData, range));
2349  CG_CHECK(cg_descriptor_write, ("TraceRange", range.c_str()));
2350  }
2351  break;
2352 
2353  case 'e':
2354  if (pane.is_structured()) {
2355  rank = cellDim;
2356  size[0] = pane.size_i() - 1;
2357  if (cellDim > 1) {
2358  size[1] = pane.size_j() - 1;
2359  size[2] = cellDim > 2 ? pane.size_k() - 1 : 1;
2360  } else {
2361  size[1] = size[2] = 1;
2362  }
2363  numItems = size[0] * size[1] * size[2];
2364  std::fill_n(rind, 2 * physDim, pane.size_of_ghost_layers());
2365  } else {
2366  rank = 1;
2367  size[0] = numItems = (*a)->size_of_items();
2368  size[1] = size[2] = 1;
2369  std::fill_n(rind, 6, 0);
2370  rind[1] = (*a)->size_of_ghost_items();
2371  }
2372 
2373  CG_CHECK(cg_sol_find_or_create, (fn, B, Z, elemName.c_str(), CellCenter,
2374  &T, errorhandle));
2375 
2376  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T, "end"));
2377  CG_CHECK(cg_gridlocation_write, (CellCenter));
2378 
2379  // Write rind data.
2380  if (writeGhost) {
2381  CG_CHECK(cg_rind_write, (rind));
2382  DEBUG_MSG("Ghost == " << rind[1]);
2383  }
2384 
2385  CG_CHECK(cg_narrays, (&offset));
2386  //++offset;
2387  // MS
2388  // If the attribute we are trying to write is already written to
2389  // the file we just skip it.
2390  //std::cout << __FILE__ << __LINE__ << std::endl;
2391  //std::cout << " Zone = " << Z << " IntegralData_t = " << T << std::endl;
2392  //std::cout << " Read offset = " << offset << std::endl;
2393  for (int ii=1; ii<=offset; ii++)
2394  {
2395  cg_array_info(ii, arrName, &dT, &tt1, &tt2);
2396  //std::cout << "offset = " << ii <<" arrayName = " << arrName << std::endl;
2397  if ((std::string(arrName)).find(std::string((*a)->name())) != std::string::npos) {
2398  //std::cout << "The request for duplicate dataitem is ignored." << std::endl;
2399  goToNextItem = true;
2400  }
2401  }
2402  if (goToNextItem)
2403  break;
2404  ++offset;
2405  //if (offset == 4)
2406  // return;
2407  // MS End
2408  //std::cout << __FILE__ << __LINE__ << std::endl;
2409 
2410  for (A=0; A<nComp; ++A) {
2411  //std::cout << __FILE__ << __LINE__ << " Component = " << A << std::endl;
2412  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T, "end"));
2413  //std::cout << __FILE__ << __LINE__ << std::endl;
2414 
2415  label = (*a)->name();
2416  //std::cout << __FILE__ << __LINE__ << std::endl;
2417 
2418  const Attribute* pa;
2419  if (nComp == 1)
2420  pa = pane.attribute((*a)->id());
2421  else {
2422  pa = pane.attribute((*a)->id() + A + 1);
2423  if (nComp == physDim) {
2424  label += (char)('X' + (char)A);
2425  } else if (nComp == physDim * physDim) {
2426  label += (char)('X' + (char)(A / physDim));
2427  label += (char)('X' + (char)(A % physDim));
2428  } else {
2429  std::ostringstream sout;
2430  sout << label << '#' << A + 1 << "of" << nComp;
2431  label = sout.str();
2432  }
2433  }
2434  //std::cout << __FILE__ << __LINE__ << std::endl;
2435 
2436 #ifdef DEBUG_DUMP_PREFIX
2437  std::ofstream fout((DEBUG_DUMP_PREFIX + std::string(material) + '.' + (*a)->name() + '_' + timeLevel + ".cgns").c_str(), std::ios_base::app);
2438 #endif // DEBUG_DUMP_PREFIX
2439  bool isNull = false;
2440  pData[A] = pa->pointer();
2441  //std::cout << __FILE__ << __LINE__ << " numItems = " << numItems << std::endl;
2442  if (numItems <= 0) {
2443  //std::cout << __FILE__ << __LINE__ << std::endl;
2444  buf.resize(COM::Attribute::get_sizeof(dType, 1), '\0');
2445  std::fill(buf.begin(), buf.end(), '\0');
2446  pData[A] = &(buf[0]);
2447  size[0] = size[1] = size[2] = 1;
2448  range = "EMPTY";
2449  } else if (pData[A] == NULL) {
2450  //std::cout << __FILE__ << __LINE__ << std::endl;
2451  buf.resize(COM::Attribute::get_sizeof(dType, std::max(numItems, 1)),
2452  '\0');
2453  std::fill(buf.begin(), buf.end(), '\0');
2454  pData[A] = &(buf[0]);
2455  isNull = true;
2456  // size[0] = size[1] = size[2] = 1;
2457  range = "NULL";
2458  } else {
2459  //std::cout << __FILE__ << __LINE__ << std::endl;
2460  int stride = pa->stride();
2461  //std::cout << __FILE__ << __LINE__ << std::endl;
2462  if (stride > 1) {
2463  //std::cout << __FILE__ << __LINE__ << std::endl;
2464  int dSize = COM::Attribute::get_sizeof(dType, 1);
2465  buf.resize(dSize * std::max(numItems, 1));
2466  //std::cout << __FILE__ << __LINE__ << std::endl;
2467  for (i=0; i<numItems; ++i) {
2468  std::memcpy(&(buf[i*dSize]),
2469  &(((const char*)pData[A])[i*stride*dSize]), dSize);
2470  }
2471  //std::cout << __FILE__ << __LINE__ << std::endl;
2472  pData[A] = &(buf[0]);
2473  }
2474 #ifdef DEBUG_DUMP_PREFIX
2475  {
2476  switch (COM2CGNS(dType)) {
2477  case Character:
2478  for (i=0; i<numItems; ++i)
2479  fout << i << " : " << (int)((const char*)(pData[A]))[i]
2480  << '\n';
2481  break;
2482  case Integer:
2483  for (i=0; i<numItems; ++i)
2484  fout << i << " : " << ((const int*)(pData[A]))[i] << '\n';
2485  break;
2486  case RealSingle:
2487  for (i=0; i<numItems; ++i)
2488  fout << i << " : " << ((const float*)(pData[A]))[i] << '\n';
2489  break;
2490  case RealDouble:
2491  for (i=0; i<numItems; ++i)
2492  fout << i << " : " << ((const double*)(pData[A]))[i]
2493  << '\n';
2494  break;
2495  default:
2496  break;
2497  }
2498  fout << "###########################################\n";
2499  }
2500 #endif // DEBUG_DUMP_PREFIX
2501 
2502  SwitchOnCOMDataType(dType,
2503  FindRange(rank, size, rind,
2504  (const COM_TT*)pData[A], range));
2505  }
2506  if (writeGhost) {
2507  DEBUG_MSG("Calling cg_array_write( name == '" << label
2508  << "', dataType == "
2509  << (COM2CGNS(dType) == RealSingle ? "float"
2510  : (COM2CGNS(dType) == RealDouble ? "double"
2511  : (COM2CGNS(dType) == Character ? "char" : "int?")))
2512  << ", rank == " << rank << ", size[] = { " << size[0]
2513  << ", " << size[1] << ", " << size[2] << " } )");
2514  CG_CHECK(cg_array_write, (label.c_str(), COM2CGNS(dType),
2515  rank, size, pData[A]));
2516  } else {
2517  DEBUG_MSG("Calling cg_array_core_write( name == '" << label
2518  << "', dataType == "
2519  << (COM2CGNS(dType) == RealSingle ? "float"
2520  : (COM2CGNS(dType) == RealDouble ? "double"
2521  : (COM2CGNS(dType) == Character ? "char" : "int?")))
2522  << ", rank == " << rank << ", rind[] = { " << rind[0]
2523  << ", " << rind[1] << ", " << rind[2] << ", " << rind[3]
2524  << ", " << rind[4] << ", " << rind[5] << " }, size[] = { "
2525  << size[0] << ", " << size[1] << ", " << size[2]
2526  << " } )");
2527  CG_CHECK(cg_array_core_write, (label.c_str(), COM2CGNS(dType),
2528  rank, rind, size, pData[A]));
2529  }
2530 
2531  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T,
2532  "DataArray_t", A + offset, "end"));
2533 
2534  DEBUG_MSG("Calling cg_descriptor_write( name == 'Range', "
2535  << "value == '" << range << "' )");
2536  CG_CHECK(cg_descriptor_write, ("Range", range.c_str()));
2537  if (!pa->unit().empty()) {
2538  cg_exponents_as_string_write(attr->unit().c_str(), errorhandle);
2539  DEBUG_MSG("Calling cg_descriptor_write( name == 'Units', "
2540  << "value == '" << pa->unit() << "' )");
2541  CG_CHECK(cg_descriptor_write, ("Units", pa->unit().c_str()));
2542  }
2543  }
2544  //std::cout << __FILE__ << __LINE__ << std::endl;
2545 
2546  // Vectors and tensors need a MagnitudeRange or TraceRange
2547  // descriptor under the first DataArray_t node.
2548  if (nComp == 3) {
2549  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T,
2550  "DataArray_t", offset, "end"));
2551  SwitchOnCOMDataType(dType,
2552  FindMagnitudeRange(physDim, rank, size, rind,
2553  (const COM_TT**)pData,
2554  range));
2555  CG_CHECK(cg_descriptor_write, ("MagnitudeRange", range.c_str()));
2556  } else if (nComp == 9) {
2557  CG_CHECK(cg_goto, (fn, B, "Zone_t", Z, "FlowSolution_t", T,
2558  "DataArray_t", offset, "end"));
2559  SwitchOnCOMDataType(dType,
2560  FindTraceRange(physDim, rank, size, rind,
2561  (const COM_TT**)pData, range));
2562  CG_CHECK(cg_descriptor_write, ("TraceRange", range.c_str()));
2563  }
2564  break;
2565  }
2566  }
2567 }
2568 
2569 
2570 
2571 
2572 
2573 
cg_base_find_or_create
static int cg_base_find_or_create(int fn, const char *name, int cellDim, int physDim, int *B, const std::string &errorhandle)
Find the named Base_t node, or create one if it doesn&#39;t exist.
Definition: Rocout_cgns.C:701
Pane::size_of_real_elements
Size size_of_real_elements() const
Get the number of real elements in the pane (excluding ghost elements).
Definition: Pane.h:158
Pane
A Pane object contains a mesh, pane attribute, and field variables.
Definition: Pane.h:43
cimg_library::exp
CImg< _cimg_Tfloat > exp(const CImg< T > &instance)
Definition: CImg.h:6016
FindTraceRange
void FindTraceRange(int physDim, int rank, const int *size, const int *rind, const T **pData, std::string &range)
Find the minimum and maximum core trace.
Definition: Rocout_cgns.C:473
FindRange
void FindRange(int rank, const int *size, const int *rind, const T *pData, std::string &range)
Find the minimum and maximum core value.
Definition: Rocout_cgns.C:342
COM_assertion
#define COM_assertion(EX)
Error checking utility similar to the assert macro of the C language.
Definition: roccom_assertion.h:62
AutoCloser::errorhandle
const std::string & errorhandle
Definition: Rocout_cgns.C:131
COM_FLOAT
Definition: roccom_basic.h:131
AutoCDer
Automatically save and restore the current working directory.
Definition: Rocout_cgns.C:139
Attribute
An Attribute object is a data member of a window.
Definition: Attribute.h:51
k
j indices k indices k
Definition: Indexing.h:6
COM_INTEGER
Definition: roccom_basic.h:133
COM_assertion_msg
#define COM_assertion_msg(EX, msg)
Definition: roccom_assertion.h:64
COM2CGNS
static DataType_t COM2CGNS(int dType)
Convert a Roccom data type to a CGNS data type.
Definition: Rocout_cgns.C:156
s
double s
Definition: blastest.C:80
FindMagnitudeRange
void FindMagnitudeRange(int physDim, int rank, const int *size, const int *rind, const T **pData, std::string &range)
Find the minimum and maximum core magnitude.
Definition: Rocout_cgns.C:418
cg_integral_find_or_create
static int cg_integral_find_or_create(const char *name, int *I, const std::string &errorhandle)
Find the named IntegralData_t node, or create one if it doesn&#39;t exist.
Definition: Rocout_cgns.C:789
AutoCDer::m_cwd
char * m_cwd
Definition: Rocout_cgns.C:147
cg_zone_find_or_create
static int cg_zone_find_or_create(int fn, int B, const char *name, const int *sizes, ZoneType_t zType, int *Z, const std::string &errorhandle)
Find the named Zone_t node, or create one if it doesn&#39;t exist.
Definition: Rocout_cgns.C:741
Rocout_cgns.h
Declaration of Rocout CGNS routines.
max
Vector_n max(const Array_n_const &v1, const Array_n_const &v2)
Definition: Vector_n.h:354
Window
A Window object contains multiple panes and multiple data attributes.
Definition: Window.h:42
cg_exponents_as_string_write
static void cg_exponents_as_string_write(std::string unit, const std::string &errorhandle)
Definition: Rocout_cgns.C:648
Pane::is_structured
bool is_structured() const
Is mesh of the pane structured?
Definition: Pane.h:110
Element_accessors.h
COM_CHAR
C/C++ Data types.
Definition: roccom_basic.h:129
Pane::size_of_ghost_nodes
Size size_of_ghost_nodes() const
Get the number of ghost nodes.
Definition: Pane.h:126
AutoCloser
Automatically close open files, datasets, etc.
Definition: Rocin.C:179
sqrt
double sqrt(double d)
Definition: double.h:73
offset
real *8 function offset(vNorm, x2, y2, z2)
Definition: PlaneNorm.f90:211
Pane::size_k
Size size_k() const
Get the number of nodes in k-dimension if the mesh is structured.
Definition: Pane.h:222
write_attr_CGNS
void write_attr_CGNS(const std::string &fname, const std::string &mfile, const COM::Attribute *attr, const char *material, const char *timelevel, int pane_id, const std::string &ghosthandle, const std::string &errorhandle, int mode)
Write the data for the given attribute to file.
Definition: Rocout_cgns.C:884
Attribute::location
Shorter_size location() const
Obtain the location of the attribute.
Definition: Attribute.h:186
ModifyExponents
void ModifyExponents(const std::string &unit, float *exponents, float dir)
Definition: Rocout_cgns.C:544
Pane::attribute
Attribute * attribute(const std::string &a)
Obtain the attribute from given name.
Definition: Pane.C:148
Attribute::unit
const std::string & unit() const
Obtain the unit of the attribute.
Definition: Attribute.h:200
Mesquite::length
double length(Vector3D *const v, int n)
Definition: includeLinks/Vector3D.hpp:400
DEBUG_MSG
#define DEBUG_MSG(x)
Definition: Rocout_cgns.C:57
Rocout.h
Rocout creates a series of Roccom windows by reading in a list of files.
Attribute::pointer
const void * pointer() const
Obtain a constant pointer to the physical address.
Definition: Attribute.h:150
A
rational * A
Definition: vinci_lass.c:67
Attribute::size_of_ghost_items
int size_of_ghost_items() const
Obtain the number of ghost items in the attribute.
Definition: Attribute.C:139
WriteElements
static int WriteElements(int fn, int B, int Z, int eType, int nElems, const int *elems, bool isStaggered, int length, int offset, const std::string &label, const std::string &errorhandle)
Build and write out an Elements_t node.
Definition: Rocout_cgns.C:201
AutoCloser::AutoCloser
AutoCloser(int fn, const std::string &eh)
Definition: Rocout_cgns.C:125
SwitchOnCOMDataType
#define SwitchOnCOMDataType(dType, funcCall)
Call a templated function based on a Roccom data type.
Definition: Rocout_cgns.C:93
Pane::size_of_ghost_layers
Size size_of_ghost_layers() const
Dimension of the pane.
Definition: Pane.h:212
sin
NT & sin
Definition: rational_rotation.h:78
cg_array_info_by_name
static int cg_array_info_by_name(const char *name, int *A, DataType_t *dType, int *rank, int *size, const std::string &errorhandle)
Return information on the named DataArray_t node.
Definition: Rocout_cgns.C:852
i
blockLoc i
Definition: read.cpp:79
Pane::size_j
Size size_j() const
Get the number of nodes in j-dimension if the mesh is structured.
Definition: Pane.h:219
x
void int int REAL * x
Definition: read.cpp:74
n
const NT & n
Definition: rational_rotation.h:72
get_sizeof
static int get_sizeof(MPI_Datatype i)
Get the size of a given MPI data type.
Definition: commpi.C:44
Attribute::stride
int stride() const
Obtain the stride of the attribute in base datatype.
Definition: Attribute.h:233
Pane::is_unstructured
bool is_unstructured() const
Is mesh of the pane unstructured?
Definition: Pane.h:103
cg_sol_find_or_create
static int cg_sol_find_or_create(int fn, int B, int Z, const char *name, GridLocation_t location, int *S, const std::string &errorhandle)
Find the named FlowSolution_t node, or create one if it doesn&#39;t exist.
Definition: Rocout_cgns.C:819
Attribute::size_of_real_items
int size_of_real_items() const
Obtain the number of real items in the attribute.
Definition: Attribute.C:167
Pane::elements
void elements(std::vector< Connectivity * > &es)
Obtain all the element connectivities of the pane.
Definition: Pane.h:268
COM_DOUBLE_PRECISION
Definition: roccom_basic.h:134
COM_INT
Definition: roccom_basic.h:130
min
Vector_n min(const Array_n_const &v1, const Array_n_const &v2)
Definition: Vector_n.h:346
ExtractExponent
float ExtractExponent(std::string &unit)
Determine and write the exponents for the basic units of measurement.
Definition: Rocout_cgns.C:527
Pane::size_of_real_nodes
Size size_of_real_nodes() const
Get the number of real nodes in the pane (excluding ghost nodes).
Definition: Pane.h:134
cg_array_core_write_internal
static int cg_array_core_write_internal(char const *name, DataType_t dType, int rank, int *rind, int const *size, T const *data)
Definition: Rocout_cgns.C:592
cwd
char cwd[FILENAME_MAX]
Definition: Rocout_cgns.C:48
Window::pane
Pane & pane(const int pane_id, bool insert=false)
Find the pane with given ID. If not found, insert a pane with given ID.
Definition: Window.C:769
j
j indices j
Definition: Indexing.h:6
Pane::size_of_nodes
Size size_of_nodes() const
Get the total number of nodes in the pane (including ghost nodes).
Definition: Pane.h:118
Pane::dimension
int dimension() const
Dimension of the pane.
Definition: Pane.h:99
COM_REAL
Definition: roccom_basic.h:133
cg_array_core_write
static int cg_array_core_write(char const *name, DataType_t dType, int rank, int *rind, int const *size, void const *data)
Definition: Rocout_cgns.C:617
GetCurrentDir
#define GetCurrentDir
Definition: Rocout_cgns.C:49
rank
static int rank
Definition: advectest.C:66
Attribute::data_type
COM_Type data_type() const
Obtain the data type of each component of the attribute.
Definition: Attribute.h:197
Attribute::id
int id() const
Obtain the id (or index) of the attribute.
Definition: Attribute.h:120
COM_CHARACTER
Fortran Data types.
Definition: roccom_basic.h:133
Pane::size_i
Size size_i() const
Get the number of nodes in i-dimension if the mesh is structured.
Definition: Pane.h:216
COM_DOUBLE
Definition: roccom_basic.h:131
Pane::attributes
void attributes(std::vector< Attribute * > &as)
Obtain all the attributes of the pane.
Definition: Pane.C:165
CG_CHECK
#define CG_CHECK(routine, args)
Perform CGNS error-checking.
Definition: Rocout_cgns.C:68
Pane::id
int id() const
Get the ID of the pane.
Definition: Pane.h:96
Attribute::size_of_items
int size_of_items() const
Obtain the number of items in the attribute.
Definition: Attribute.C:111
Attribute::size_of_components
int size_of_components() const
Obtain the number of components in the attribute.
Definition: Attribute.h:203