ccmReaderSolution.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2016-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*----------------------------------------------------------------------------*/
27
28#include "ccmReader.H"
29#include "ccmInternal.H" // include last to avoid any strange interactions
30
31
32// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33
34// Get information about available fields.
35// - assume that all fields are available for all solution intervals
36// eg,
37// /FieldSets/FieldSet-1/Phase-1/Fields/Velocity
38// /FieldSets/FieldSet-1/Phase-1/Fields/Pressure
39// ...
40void Foam::ccm::reader::determineFieldInfo
41(
42 const ccmID& fieldSetNode,
43 fieldTable& table
44)
45{
46 char fullName[kCCMIOMaxStringLength + 1];
47 char shortName[kCCMIOProstarShortNameLength+1];
48
49 // Loop through all phases
50 ccmID phaseNode;
51 int phaseI = 0;
52
53 while
54 (
55 CCMIONextEntity
56 (
57 nullptr,
58 fieldSetNode,
59 kCCMIOFieldPhase,
60 &phaseI,
61 &phaseNode
62 )
63 == kCCMIONoErr
64 )
65 {
66 CCMIODimensionality dims;
67
68 // Examine each field in this fieldSet
69 ccmID fieldNode;
70 int fieldI = 0;
71
72 // Get full/short names and dimension (scalar/vector/tensor).
73 // Use short name as unique identifier.
74
75 while
76 (
77 CCMIONextEntity
78 (
79 nullptr,
80 phaseNode,
81 kCCMIOField,
82 &fieldI,
83 &fieldNode
84 )
85 == kCCMIONoErr
86
87 && CCMIOReadField
88 (
89 nullptr,
90 fieldNode,
91 fullName,
92 shortName,
93 &dims,
94 nullptr
95 )
96 == kCCMIONoErr
97 )
98 {
99 // The shortName *should* be a word, but some fields seem
100 // to have blanks!
101 {
102 char *ptr = shortName;
103 while (*ptr)
104 {
105 if (!word::valid(*ptr))
106 {
107 *ptr = '_';
108 }
109 ptr++;
110 }
111 }
112
113 word fieldName(shortName);
114
115 if (dims == kCCMIOScalar)
116 {
117 // Add to table as required
118 if (!table.found(fieldName))
119 {
120 fieldEntry entry(fieldName, fullName);
121
122 entry.units
123 (
124 ccmReadOptstr("Units", fieldNode)
125 );
126
127 table.append(entry);
128 }
129
130 fieldEntry& entry = table.find(fieldName)();
131
132 // Obtain sizes of data field
133 // - only process cell/face data
134 ccmID dataNode;
135 int dataI = 0;
136
137 CCMIODataLocation dataLocation;
138 CCMIOIndex maxId;
139
140 while
141 (
142 CCMIONextEntity
143 (
144 nullptr,
145 fieldNode,
146 kCCMIOFieldData,
147 &dataI,
148 &dataNode
149 )
150 == kCCMIONoErr
151
152 && CCMIOEntitySize
153 (
154 nullptr,
155 dataNode,
156 nullptr,
157 &maxId
158 )
159 == kCCMIONoErr
160
161 && CCMIOReadFieldDatad
162 (
163 nullptr,
164 dataNode,
165 nullptr,
166 &dataLocation,
167 nullptr,
168 kCCMIOStart,
169 kCCMIOEnd
170 )
171 == kCCMIONoErr
172 )
173 {
174 if (dataLocation == kCCMIOCell)
175 {
176 entry.maxCellId(maxId);
177 }
178 else if (dataLocation == kCCMIOFace)
179 {
180 entry.maxFaceId(maxId);
181 }
182 }
183 }
184 }
185 }
186}
187
188
189// Get information about all available solutions.
190// - it is sufficient to examine the first processor
191//
192// The "States":
193// Restart_1/Processor-1
194// -OR-
195// Transient_1/Processor-1
196// ...
197// Transient_N/Processor-1
198//
199// point to the "FieldSets":
200// FieldSet-1/RestartInfo
201//
202bool Foam::ccm::reader::detectSolution()
203{
204 // call once
205 if (solutionStatus_ != UNKNOWN)
206 {
207 return (solutionStatus_ == OKAY || solutionStatus_ == READ);
208 }
209
210 // loop through all States
211 ccmID stateNode;
212 int stateI = 0;
213 while
214 (
215 CCMIONextEntity
216 (
217 nullptr,
218 (globalState_->root),
219 kCCMIOState,
220 &stateI,
221 &stateNode
222 )
223 == kCCMIONoErr
224 )
225 {
226 // Loop through all processors/solutions
227 ccmID processorNode;
228 ccmID solutionNode;
229 int procI = 0;
230
231 while
232 (
233 CCMIONextEntity
234 (
235 nullptr,
236 stateNode,
237 kCCMIOProcessor,
238 &procI,
239 &processorNode
240 )
241 == kCCMIONoErr
242
243 && CCMIOReadProcessor
244 (
245 nullptr,
246 processorNode,
247 nullptr, // Ignore verticesNode
248 nullptr, // Ignore topologyNode
249 nullptr, // Ignore initialField
250 &solutionNode
251 )
252 == kCCMIONoErr
253 )
254 {
255 // Restrict to solutions with RestartInfo on the first cpu
256 // (there is normally only one set of restart data)
257 ccmID restartNode;
258 int restartI = 0;
259
260 char solutionName[kCCMIOMaxStringLength + 1];
261 int iteration = 0;
262 float timeValue = 0;
263
264 if
265 (
266 CCMIONextEntity
267 (
268 nullptr,
269 solutionNode,
270 kCCMIORestart,
271 &restartI,
272 &restartNode
273 )
274 == kCCMIONoErr
275
276 && CCMIOEntityName
277 (
278 nullptr,
279 stateNode,
280 solutionName
281 )
282 == kCCMIONoErr
283
284 && CCMIOReadRestartInfo
285 (
286 nullptr,
287 restartNode,
288 nullptr, // Ignore solverName
289 &iteration,
290 &timeValue,
291 nullptr, // Ignore timeUnits
292 nullptr // Ignore startAngle
293 )
294 == kCCMIONoErr
295 )
296 {
297 solutionTable_.append
298 (
299 solutionEntry(solutionName, iteration, timeValue)
300 );
301 }
302
303 // Determine field information
304 determineFieldInfo(solutionNode, fieldTable_);
305
306 // check Lagrangian data
307 ccmID lagrangianNode;
308 ccmID lagrangianSolutions;
309 int lagrangianI = 0;
310
311 if
312 (
313 CCMIONextEntity
314 (
315 nullptr,
316 processorNode,
317 kCCMIOLagrangianData,
318 &lagrangianI,
319 &lagrangianNode
320 )
321 == kCCMIONoErr
322
323 && CCMIOReadLagrangianData
324 (
325 nullptr,
326 lagrangianNode,
327 nullptr,
328 &lagrangianSolutions
329 )
330 == kCCMIONoErr
331 )
332 {
333 determineFieldInfo(lagrangianSolutions, lagrangianTable_);
334 }
335 }
336 }
337
338 if (solutionTable_.size() && fieldTable_.size())
339 {
340 solutionStatus_ = OKAY;
341 }
342 else
343 {
344 solutionStatus_ = BAD;
345 }
346
347 return (solutionStatus_ == OKAY || solutionStatus_ == READ);
348}
349
350#define SOLID_STRESS_HACK
351
352//
353// Read solution and field combination
354//
357(
358 const word& solutionName,
359 const word& fieldName,
360 const bool wallData
361)
362{
363 // get State by name
364 ccmID stateNode;
365 ccmID solutionNode;
366
367 if
368 (
369 CCMIOGetState
370 (
371 nullptr,
372 (globalState_->root),
373 solutionName.c_str(),
374 nullptr,
375 &stateNode
376 )
377 != kCCMIONoErr
378
379 || !fieldTable_.found(fieldName)
380 )
381 {
382 return tmp<scalarField>::New();
383 }
384
385 CCMIODataLocation requestedLocation = kCCMIOCell;
386
387 fieldEntry& entry = fieldTable_.find(fieldName)();
388
389 label maxId = entry.maxCellId();
390
391 if (wallData)
392 {
393 maxId = entry.maxFaceId();
394 requestedLocation = kCCMIOFace;
395 }
396
397 // we can skip empty fields immediately
398 if (!maxId)
399 {
400 return tmp<scalarField>::New();
401 }
402
403 char shortName[kCCMIOProstarShortNameLength+1];
404
405 List<label> mapData;
406 List<scalar> rawData;
407
408 tmp<scalarField> tscalarData
409 (
410 new scalarField(maxId + 1, option().undefScalar())
411 );
412 scalarField& scalarData = tscalarData.ref();
413
414
415 CCMIODimensionality dims;
416
417 // Loop across all processors
418 ccmID processorNode;
419 int procI = 0;
420
421 while
422 (
423 CCMIONextEntity
424 (
425 nullptr,
426 stateNode,
427 kCCMIOProcessor,
428 &procI,
429 &processorNode
430 )
431 == kCCMIONoErr
432
433 && CCMIOReadProcessor
434 (
435 nullptr,
436 processorNode,
437 nullptr, // Ignore verticesNode
438 nullptr, // Ignore topologyNode
439 nullptr, // Ignore initialField
440 &solutionNode
441 )
442 == kCCMIONoErr
443 )
444 {
445 // loop through all phases
446 int phaseI = 0;
447 ccmID phaseNode;
448
449 while
450 (
451 CCMIONextEntity
452 (
453 nullptr,
454 solutionNode,
455 kCCMIOFieldPhase,
456 &phaseI,
457 &phaseNode
458 )
459 == kCCMIONoErr
460 )
461 {
462 // Get Field that matches the name
463 ccmID fieldNode;
464 int fieldI = 0;
465
466 while
467 (
468 CCMIONextEntity
469 (
470 nullptr,
471 phaseNode,
472 kCCMIOField,
473 &fieldI,
474 &fieldNode
475 )
476 == kCCMIONoErr
477 )
478 {
479 // Get full/short names and dimension (scalar/vector/tensor)
480 // use short name as unique identifier
481 CCMIOReadField
482 (
483 &(globalState_->error),
484 fieldNode,
485 nullptr,
486 shortName,
487 &dims,
488 nullptr
489 );
490 assertNoError
491 (
492 "reading post data field: "
493 + string(shortName)
494 );
495
496 if (fieldName == shortName && dims == kCCMIOScalar)
497 {
498 // Obtain sizes of data field
499 ccmID dataNode;
500 ccmID mapId;
501 int dataI = 0;
502 while
503 (
504 CCMIONextEntity
505 (
506 nullptr,
507 fieldNode,
508 kCCMIOFieldData,
509 &dataI,
510 &dataNode
511 )
512 == kCCMIONoErr
513 )
514 {
515 CCMIODataLocation dataLocation;
516 CCMIOSize n;
517
518 // Only process cell/face data
519 if
520 (
521 CCMIOEntitySize
522 (
523 nullptr,
524 dataNode,
525 &n,
526 nullptr
527 )
528 == kCCMIONoErr
529
530 && CCMIOReadFieldDatad
531 (
532 nullptr,
533 dataNode,
534 &mapId,
535 &dataLocation,
536 nullptr,
537 kCCMIOStart,
538 kCCMIOEnd
539 )
540 == kCCMIONoErr
541
542 && dataLocation == requestedLocation
543 )
544 {
545#ifdef SOLID_STRESS_HACK
546 bool okayCombination = true;
547
548 CCMIOSize len;
549 if
550 (
551 CCMIOEntityDescription
552 (
553 nullptr,
554 dataNode,
555 &len,
556 nullptr
557 )
558 == kCCMIONoErr
559 )
560 {
561 char* dataLabel = new char[len + 1];
562
563 if
564 (
565 CCMIOEntityDescription
566 (
567 nullptr,
568 dataNode,
569 &len,
570 dataLabel
571 )
572 == kCCMIONoErr
573 )
574 {
575 if
576 (
577 (
578 strstr(fieldName.c_str(), "SIG")
579 || strstr(fieldName.c_str(), "EPS")
580 )
581 && strstr(dataLabel, "So") == nullptr
582 )
583 {
584 okayCombination = false;
585 // skip non-solid
586 }
587 }
588
589 delete[] dataLabel;
590 }
591
592 if (!okayCombination)
593 {
594 continue;
595 }
596#endif
597
598 mapData.setSize(n);
599 rawData.setSize(n);
600
601 readMap
602 (
603 mapId,
604 mapData
605 );
606
607 CCMIOReadFieldDatad
608 (
609 &(globalState_->error),
610 dataNode,
611 nullptr,
612 nullptr,
613 rawData.data(),
614 kCCMIOStart,
615 kCCMIOEnd
616 );
617 assertNoError
618 (
619 "reading post data field: "
620 + string(shortName)
621 );
622
623 // transcribe to output list
624 forAll(mapData, i)
625 {
626 const label cellId = mapData[i];
627 scalarData[cellId] = rawData[i];
628 }
629 }
630 }
631 }
632 }
633 }
634 }
635
636 // Overwrite possible junk in cell 0 (often used as /dev/null)
637 scalarData[0] = option().undefScalar();
638
639 return tscalarData;
640}
641
642
643// ************************************************************************* //
label n
Internal bits for wrapping libccmio - do not use directly.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: List.H:77
void setSize(const label n)
Alias for resize()
Definition: List.H:218
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
T * data() noexcept
Return pointer to the underlying array serving as data storage.
Definition: UListI.H:237
A ccm field entry with short name, name, maxId and type.
bool valid() const
True if all internal ids are non-negative.
tmp< scalarField > readField(const word &solutionName, const word &fieldName, const bool wallData=false)
Read solution and field combination.
A keyword and a list of tokens is an 'entry'.
Definition: entry.H:70
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
A class for handling words, derived from Foam::string.
Definition: word.H:68
label cellId
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333