"Fossies" - the Fresh Open Source Software Archive  

Source code changes of the file "Kernel/Pos_Format.cpp" between
getdp-3.4.0-source.tgz and getdp-3.5.0-source.tgz

About: GetDP is a general finite element solver using mixed elements to discretize de Rham-type complexes in one, two and three dimensions.

Pos_Format.cpp  (getdp-3.4.0-source.tgz):Pos_Format.cpp  (getdp-3.5.0-source.tgz)
// GetDP - Copyright (C) 1997-2021 P. Dular and C. Geuzaine, University of Liege // GetDP - Copyright (C) 1997-2022 P. Dular and C. Geuzaine, University of Liege
// //
// See the LICENSE.txt file for license information. Please report all // See the LICENSE.txt file for license information. Please report all
// issues on https://gitlab.onelab.info/getdp/getdp/issues. // issues on https://gitlab.onelab.info/getdp/getdp/issues.
#include <sstream> #include <sstream>
#include <map> #include <map>
#include <vector> #include <vector>
#include <list> #include <list>
#include <string.h> #include <string.h>
#include <math.h> #include <math.h>
skipping to change at line 38 skipping to change at line 38
#include "MallocUtils.h" #include "MallocUtils.h"
#include "Message.h" #include "Message.h"
#include "kissfft.hh" #include "kissfft.hh"
#if defined(HAVE_GMSH) #if defined(HAVE_GMSH)
#include <gmsh/PView.h> #include <gmsh/PView.h>
#include <gmsh/PViewDataList.h> #include <gmsh/PViewDataList.h>
#include <gmsh/PViewDataGModel.h> #include <gmsh/PViewDataGModel.h>
#endif #endif
#define TWO_PI 6.2831853071795865 #define TWO_PI 6.2831853071795865
#define NBR_MAX_ISO 200 #define NBR_MAX_ISO 200
#define SQU(a) ((a)*(a)) #define SQU(a) ((a) * (a))
extern struct Problem Problem_S ; extern struct Problem Problem_S;
extern struct CurrentData Current ; extern struct CurrentData Current;
extern int Flag_BIN, Flag_GMSH_VERSION; extern int Flag_BIN, Flag_GMSH_VERSION;
extern FILE *PostStream ; extern FILE *PostStream;
static List_T *PostElement_L = NULL ; static List_T *PostElement_L = NULL;
static List_T *TimeValue_L = NULL ; static List_T *TimeValue_L = NULL;
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Gmsh formats */ /* Gmsh formats */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
// global static lists for new-style Gmsh output (cannot be saved incrementally // global static lists for new-style Gmsh output (cannot be saved incrementally
// for each element) // for each element)
static int Gmsh_StartNewView = 0 ; static int Gmsh_StartNewView = 0;
static int NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT, NbTT; static int NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT, NbTT;
static int NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH; static int NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH;
static int NbSI, NbVI, NbTI, NbSY, NbVY, NbTY; static int NbSI, NbVI, NbTI, NbSY, NbVY, NbTY;
static int NbT2; static int NbT2;
static std::vector<double> SP, VP, TP, SL, VL, TL, ST, VT, TT, SQ, VQ, TQ; static std::vector<double> SP, VP, TP, SL, VL, TL, ST, VT, TT, SQ, VQ, TQ;
static std::vector<double> SS, VS, TS1 /* for petsc */, SH, VH, TH; static std::vector<double> SS, VS, TS1 /* for petsc */, SH, VH, TH;
static std::vector<double> SI, VI, TI, SY, VY, TY, T2D; static std::vector<double> SI, VI, TI, SY, VY, TY, T2D;
static std::vector<char> T2C; static std::vector<char> T2C;
static char CurrentName[256] = ""; static char CurrentName[256] = "";
static int CurrentPartitionNumber = 0; static int CurrentPartitionNumber = 0;
static void Gmsh_ResetStaticLists() static void Gmsh_ResetStaticLists()
{ {
NbSP = NbVP = NbTP = NbSL = NbVL = NbTL = 0; NbSP = NbVP = NbTP = NbSL = NbVL = NbTL = 0;
NbST = NbVT = NbTT = NbSQ = NbVQ = NbTQ = 0; NbST = NbVT = NbTT = NbSQ = NbVQ = NbTQ = 0;
NbSS = NbVS = NbTS = NbSH = NbVH = NbTH = 0; NbSS = NbVS = NbTS = NbSH = NbVH = NbTH = 0;
NbSI = NbVI = NbTI = NbSY = NbVY = NbTY = 0; NbSI = NbVI = NbTI = NbSY = NbVY = NbTY = 0;
NbT2 = 0; NbT2 = 0;
SP.clear(); VP.clear(); TP.clear(); SL.clear(); VL.clear(); TL.clear(); SP.clear();
ST.clear(); VT.clear(); TT.clear(); SQ.clear(); VQ.clear(); TQ.clear(); VP.clear();
SS.clear(); VS.clear(); TS1.clear(); SH.clear(); VH.clear(); TH.clear(); TP.clear();
SI.clear(); VI.clear(); TI.clear(); SY.clear(); VY.clear(); TY.clear(); SL.clear();
T2D.clear(); T2C.clear(); VL.clear();
if(!TimeValue_L) TimeValue_L = List_Create(100,1000000,sizeof(double)); TL.clear();
else List_Reset(TimeValue_L); ST.clear();
VT.clear();
TT.clear();
SQ.clear();
VQ.clear();
TQ.clear();
SS.clear();
VS.clear();
TS1.clear();
SH.clear();
VH.clear();
TH.clear();
SI.clear();
VI.clear();
TI.clear();
SY.clear();
VY.clear();
TY.clear();
T2D.clear();
T2C.clear();
if(!TimeValue_L)
TimeValue_L = List_Create(100, 1000000, sizeof(double));
else
List_Reset(TimeValue_L);
} }
static void Gmsh_StringStart(int Format, double x, double y, double style) static void Gmsh_StringStart(int Format, double x, double y, double style)
{ {
if(Flag_BIN){ /* bricolage: should use Format instead */ if(Flag_BIN) { /* bricolage: should use Format instead */
T2D.push_back(x); T2D.push_back(x);
T2D.push_back(y); T2D.push_back(y);
T2D.push_back(style); T2D.push_back(style);
T2D.push_back(T2C.size()); T2D.push_back(T2C.size());
NbT2++; NbT2++;
} }
else if(PostStream){ else if(PostStream) {
fprintf(PostStream, "T2(%g,%g,%g){", x, y, style); fprintf(PostStream, "T2(%g,%g,%g){", x, y, style);
} }
} }
static void Gmsh_StringAdd(int Format, int first, char *text) static void Gmsh_StringAdd(int Format, int first, char *text)
{ {
int i; int i;
if(Flag_BIN){ /* bricolage: should use Format instead */ if(Flag_BIN) { /* bricolage: should use Format instead */
for(i = 0; i < (int)strlen(text)+1; i++) for(i = 0; i < (int)strlen(text) + 1; i++) T2C.push_back(text[i]);
T2C.push_back(text[i]); }
} else if(PostStream) {
else if(PostStream){ if(!first) fprintf(PostStream, ",");
if(!first)
fprintf(PostStream, ",");
fprintf(PostStream, "\"%s\"", text); fprintf(PostStream, "\"%s\"", text);
} }
} }
static void Gmsh_StringEnd(int Format) static void Gmsh_StringEnd(int Format)
{ {
if(Flag_BIN){ /* bricolage: should use Format instead */ if(Flag_BIN) { /* bricolage: should use Format instead */
} }
else if(PostStream){ else if(PostStream) {
fprintf(PostStream, "};\n") ; fprintf(PostStream, "};\n");
} }
} }
static void Gmsh_PrintElementNodeData(struct PostSubOperation *PSO_P, int numTim static void Gmsh_PrintElementNodeData(struct PostSubOperation *PSO_P,
eStep, int numTimeStep, int numComp, int Nb[8],
int numComp, int Nb[8], std::vector<double std::vector<double> *L[8])
> *L[8])
{ {
if(!PostStream) return; if(!PostStream) return;
int N = 0; int N = 0;
for(int i = 0; i < 8; i++) N += Nb[i]; for(int i = 0; i < 8; i++) N += Nb[i];
if(!N) return; if(!N) return;
int step = 0; int step = 0;
for (int ts = 0; ts < numTimeStep; ts++) { for(int ts = 0; ts < numTimeStep; ts++) {
Pos_InitAllSolutions(PSO_P->TimeStep_L, ts); Pos_InitAllSolutions(PSO_P->TimeStep_L, ts);
for(int har = 0; har < Current.NbrHar; har++){ for(int har = 0; har < Current.NbrHar; har++) {
fprintf(PostStream, "$ElementNodeData\n"); fprintf(PostStream, "$ElementNodeData\n");
fprintf(PostStream, "1\n"); fprintf(PostStream, "1\n");
fprintf(PostStream, "\"%s\"\n", CurrentName); fprintf(PostStream, "\"%s\"\n", CurrentName);
fprintf(PostStream, "1\n"); fprintf(PostStream, "1\n");
fprintf(PostStream, "%.16g\n", Current.Time); fprintf(PostStream, "%.16g\n", Current.Time);
fprintf(PostStream, "4\n"); fprintf(PostStream, "4\n");
fprintf(PostStream, "%d\n", (PSO_P->OverrideTimeStepValue >= 0) ? fprintf(PostStream, "%d\n",
PSO_P->OverrideTimeStepValue : (Current.NbrHar > 1 ? step : (PSO_P->OverrideTimeStepValue >= 0) ?
(int)Current.TimeStep)); PSO_P->OverrideTimeStepValue :
(Current.NbrHar > 1 ? step : (int)Current.TimeStep));
fprintf(PostStream, "%d\n", numComp); fprintf(PostStream, "%d\n", numComp);
fprintf(PostStream, "%d\n", N); fprintf(PostStream, "%d\n", N);
fprintf(PostStream, "%d\n", CurrentPartitionNumber); fprintf(PostStream, "%d\n", CurrentPartitionNumber);
for(int i = 0; i < 8; i++){ for(int i = 0; i < 8; i++) {
if(!Nb[i]) continue; if(!Nb[i]) continue;
int stride = (*L[i]).size() / Nb[i]; int stride = (*L[i]).size() / Nb[i];
for(unsigned int j = 0; j < (*L[i]).size(); j += stride){ for(unsigned int j = 0; j < (*L[i]).size(); j += stride) {
double *tmp = &(*L[i])[j]; double *tmp = &(*L[i])[j];
int num = (int)tmp[0]; int num = (int)tmp[0];
int mult = (stride - 1) / numTimeStep / Current.NbrHar / numComp; int mult = (stride - 1) / numTimeStep / Current.NbrHar / numComp;
if(Flag_BIN){ if(Flag_BIN) {
fwrite(&num, sizeof(int), 1, PostStream); fwrite(&num, sizeof(int), 1, PostStream);
fwrite(&mult, sizeof(int), 1, PostStream); fwrite(&mult, sizeof(int), 1, PostStream);
fwrite(&tmp[1 + step * mult * numComp], sizeof(double), mult * numCo fwrite(&tmp[1 + step * mult * numComp], sizeof(double),
mp, mult * numComp, PostStream);
PostStream);
} }
else{ else {
fprintf(PostStream, "%d %d", num, mult); fprintf(PostStream, "%d %d", num, mult);
for(int k = 0; k < mult * numComp; k++) for(int k = 0; k < mult * numComp; k++)
fprintf(PostStream, " %.16g", tmp[1 + step * mult * numComp + k]); fprintf(PostStream, " %.16g", tmp[1 + step * mult * numComp + k]);
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
} }
} }
} }
fprintf(PostStream, "$EndElementNodeData\n"); fprintf(PostStream, "$EndElementNodeData\n");
step++; step++;
} }
} }
} }
static void GmshParsed_PrintElement(double Time, int TimeStep, int NbTimeStep, i static void GmshParsed_PrintElement(double Time, int TimeStep, int NbTimeStep,
nt NbHarmonic, int NbHarmonic, int HarmonicToTime,
int HarmonicToTime, int Type, int NbrNodes, int Type, int NbrNodes, double *x,
double *x, double *y, double *z, struct Valu double *y, double *z, struct Value *Value)
e *Value) {
{ int i, j, k, jj;
int i,j,k,jj ; double TimeMH;
double TimeMH ; struct Value TmpValue;
struct Value TmpValue ; int symIndex[9] = {0, 1, 2, 1, 3, 4, 2, 4, 5};
int symIndex[9] = {0, 1, 2, 1, 3, 4, 2, 4, 5} ; int diagIndex[9] = {0, -1, -1, -1, 1, -1, -1, -1, 2};
int diagIndex[9] = {0, -1, -1, -1, 1, -1, -1, -1, 2} ;
if(Gmsh_StartNewView){ if(Gmsh_StartNewView) {
Gmsh_StartNewView = 0 ; Gmsh_StartNewView = 0;
Gmsh_ResetStaticLists(); Gmsh_ResetStaticLists();
} }
if (HarmonicToTime == 1){ if(HarmonicToTime == 1) {
if(NbHarmonic == 2 && NbTimeStep == 1){ // classical complex case if(NbHarmonic == 2 && NbTimeStep == 1) { // classical complex case
double zero = 0., one = 1.; double zero = 0., one = 1.;
List_Put(TimeValue_L, 0, &zero); List_Put(TimeValue_L, 0, &zero);
List_Put(TimeValue_L, 1, &one); List_Put(TimeValue_L, 1, &one);
} }
else{ else {
for(k = 0 ; k < NbHarmonic ; k++) for(k = 0; k < NbHarmonic; k++)
List_Put(TimeValue_L, NbHarmonic*TimeStep+k, &Time); List_Put(TimeValue_L, NbHarmonic * TimeStep + k, &Time);
} }
} }
else else
for(k = 0 ; k < HarmonicToTime ; k++) for(k = 0; k < HarmonicToTime; k++)
List_Put(TimeValue_L, HarmonicToTime*TimeStep+k, &Time); List_Put(TimeValue_L, HarmonicToTime * TimeStep + k, &Time);
if(!PostStream) return; if(!PostStream) return;
switch (Value[0].Type) { switch(Value[0].Type) {
case SCALAR:
case SCALAR : if(TimeStep == 0) {
switch(Type) {
if(TimeStep == 0){ case POINT_ELEMENT: fprintf(PostStream, "SP("); break;
switch(Type){ case LINE: fprintf(PostStream, "SL("); break;
case POINT_ELEMENT : fprintf(PostStream, "SP("); break; case TRIANGLE: fprintf(PostStream, "ST("); break;
case LINE : fprintf(PostStream, "SL("); break; case QUADRANGLE: fprintf(PostStream, "SQ("); break;
case TRIANGLE : fprintf(PostStream, "ST("); break; case TETRAHEDRON: fprintf(PostStream, "SS("); break;
case QUADRANGLE : fprintf(PostStream, "SQ("); break; case HEXAHEDRON: fprintf(PostStream, "SH("); break;
case TETRAHEDRON : fprintf(PostStream, "SS("); break; case PRISM: fprintf(PostStream, "SI("); break;
case HEXAHEDRON : fprintf(PostStream, "SH("); break; case PYRAMID: fprintf(PostStream, "SY("); break;
case PRISM : fprintf(PostStream, "SI("); break; }
case PYRAMID : fprintf(PostStream, "SY("); break; for(i = 0; i < NbrNodes; i++) {
} if(i) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ fprintf(PostStream, "%.16g,%.16g,%.16g", x[i], y[i], z[i]);
if(i) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g,%.16g,%.16g", x[i], y[i], z[i]);
} }
fprintf(PostStream, "){"); fprintf(PostStream, "){");
} }
if (HarmonicToTime == 1) { if(HarmonicToTime == 1) {
for(k = 0 ; k < NbHarmonic ; k++) { for(k = 0; k < NbHarmonic; k++) {
if(k || TimeStep) fprintf(PostStream, ","); if(k || TimeStep) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
if(i) fprintf(PostStream, ","); if(i) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM*k]); fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM * k]);
} }
} }
} }
else { else {
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
if(k || TimeStep) fprintf(PostStream, ","); if(k || TimeStep) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
if(i) fprintf(PostStream, ","); if(i) fprintf(PostStream, ",");
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
fprintf(PostStream, "%.16g", TmpValue.Val[0]); fprintf(PostStream, "%.16g", TmpValue.Val[0]);
} }
} }
} }
if(TimeStep == NbTimeStep-1){ if(TimeStep == NbTimeStep - 1) { fprintf(PostStream, "};\n"); }
fprintf(PostStream, "};\n") ; break;
}
break ; case VECTOR:
case VECTOR : if(TimeStep == 0) {
switch(Type) {
if(TimeStep == 0){ case POINT_ELEMENT: fprintf(PostStream, "VP("); break;
switch(Type){ case LINE: fprintf(PostStream, "VL("); break;
case POINT_ELEMENT : fprintf(PostStream, "VP("); break; case TRIANGLE: fprintf(PostStream, "VT("); break;
case LINE : fprintf(PostStream, "VL("); break; case QUADRANGLE: fprintf(PostStream, "VQ("); break;
case TRIANGLE : fprintf(PostStream, "VT("); break; case TETRAHEDRON: fprintf(PostStream, "VS("); break;
case QUADRANGLE : fprintf(PostStream, "VQ("); break; case HEXAHEDRON: fprintf(PostStream, "VH("); break;
case TETRAHEDRON : fprintf(PostStream, "VS("); break; case PRISM: fprintf(PostStream, "VI("); break;
case HEXAHEDRON : fprintf(PostStream, "VH("); break; case PYRAMID: fprintf(PostStream, "VY("); break;
case PRISM : fprintf(PostStream, "VI("); break; }
case PYRAMID : fprintf(PostStream, "VY("); break; for(i = 0; i < NbrNodes; i++) {
} if(i) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ fprintf(PostStream, "%.16g,%.16g,%.16g", x[i], y[i], z[i]);
if(i) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g,%.16g,%.16g", x[i], y[i], z[i]);
} }
fprintf(PostStream, "){"); fprintf(PostStream, "){");
} }
if (HarmonicToTime == 1) { if(HarmonicToTime == 1) {
for(k = 0 ; k < NbHarmonic ; k++) { for(k = 0; k < NbHarmonic; k++) {
if(k || TimeStep) fprintf(PostStream, ","); if(k || TimeStep) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
if(i) fprintf(PostStream, ","); if(i) fprintf(PostStream, ",");
for(j = 0 ; j < 3 ; j++){ for(j = 0; j < 3; j++) {
if(j) fprintf(PostStream, ","); if(j) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM*k+j]); fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM * k + j]);
} }
} }
} }
} }
else { else {
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
if(k || TimeStep) fprintf(PostStream, ","); if(k || TimeStep) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
if(i) fprintf(PostStream, ","); if(i) fprintf(PostStream, ",");
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
for(j = 0 ; j < 3 ; j++){ for(j = 0; j < 3; j++) {
if(j) fprintf(PostStream, ","); if(j) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g", TmpValue.Val[j]); fprintf(PostStream, "%.16g", TmpValue.Val[j]);
} }
} }
} }
} }
if(TimeStep == NbTimeStep-1){ if(TimeStep == NbTimeStep - 1) { fprintf(PostStream, "};\n"); }
fprintf(PostStream, "};\n") ; break;
}
break ; case TENSOR_DIAG:
case TENSOR_SYM:
case TENSOR_DIAG : case TENSOR:
case TENSOR_SYM :
case TENSOR : if(TimeStep == 0) {
switch(Type) {
if(TimeStep == 0){ case POINT_ELEMENT: fprintf(PostStream, "TP("); break;
switch(Type){ case LINE: fprintf(PostStream, "TL("); break;
case POINT_ELEMENT : fprintf(PostStream, "TP("); break; case TRIANGLE: fprintf(PostStream, "TT("); break;
case LINE : fprintf(PostStream, "TL("); break; case QUADRANGLE: fprintf(PostStream, "TQ("); break;
case TRIANGLE : fprintf(PostStream, "TT("); break; case TETRAHEDRON: fprintf(PostStream, "TS("); break;
case QUADRANGLE : fprintf(PostStream, "TQ("); break; case HEXAHEDRON: fprintf(PostStream, "TH("); break;
case TETRAHEDRON : fprintf(PostStream, "TS("); break; case PRISM: fprintf(PostStream, "TI("); break;
case HEXAHEDRON : fprintf(PostStream, "TH("); break; case PYRAMID: fprintf(PostStream, "TY("); break;
case PRISM : fprintf(PostStream, "TI("); break; }
case PYRAMID : fprintf(PostStream, "TY("); break; for(i = 0; i < NbrNodes; i++) {
} if(i) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ fprintf(PostStream, "%.16g,%.16g,%.16g", x[i], y[i], z[i]);
if(i) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g,%.16g,%.16g", x[i], y[i], z[i]);
} }
fprintf(PostStream, "){"); fprintf(PostStream, "){");
} }
if (HarmonicToTime == 1) { if(HarmonicToTime == 1) {
for(k = 0 ; k < NbHarmonic ; k++) { for(k = 0; k < NbHarmonic; k++) {
if(k || TimeStep) fprintf(PostStream, ","); if(k || TimeStep) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
if(i) fprintf(PostStream, ","); if(i) fprintf(PostStream, ",");
for(j = 0 ; j < 9 ; j++){ for(j = 0; j < 9; j++) {
if(j) fprintf(PostStream, ","); if(j) fprintf(PostStream, ",");
if(Value[0].Type != TENSOR_DIAG) { if(Value[0].Type != TENSOR_DIAG) {
if(Value[0].Type == TENSOR_SYM) jj = symIndex[j]; if(Value[0].Type == TENSOR_SYM)
else jj = j; jj = symIndex[j];
fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM*k+jj]); else
jj = j;
fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM * k + jj]);
} }
else { else {
jj = diagIndex[j]; jj = diagIndex[j];
if(jj == -1) fprintf(PostStream, "%.16g", 0.); if(jj == -1)
else fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM*k+jj]); fprintf(PostStream, "%.16g", 0.);
else
fprintf(PostStream, "%.16g", Value[i].Val[MAX_DIM * k + jj]);
} }
} }
} }
} }
} }
else { else {
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
if(k || TimeStep) fprintf(PostStream, ","); if(k || TimeStep) fprintf(PostStream, ",");
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
if(i) fprintf(PostStream, ","); if(i) fprintf(PostStream, ",");
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
for(j = 0 ; j < 9 ; j++){ for(j = 0; j < 9; j++) {
if(j) fprintf(PostStream, ","); if(j) fprintf(PostStream, ",");
if(Value[0].Type != TENSOR_DIAG) { if(Value[0].Type != TENSOR_DIAG) {
if(Value[0].Type == TENSOR_SYM) jj = symIndex[j]; if(Value[0].Type == TENSOR_SYM)
else jj = j; jj = symIndex[j];
else
jj = j;
jj = symIndex[j]; jj = symIndex[j];
fprintf(PostStream, "%.16g", TmpValue.Val[jj]); fprintf(PostStream, "%.16g", TmpValue.Val[jj]);
} }
else { else {
jj = diagIndex[j]; jj = diagIndex[j];
if(jj == -1) fprintf(PostStream, "%.16g", 0.); if(jj == -1)
else fprintf(PostStream, "%.16g", TmpValue.Val[jj]); fprintf(PostStream, "%.16g", 0.);
else
fprintf(PostStream, "%.16g", TmpValue.Val[jj]);
} }
} }
} }
} }
} }
if(TimeStep == NbTimeStep-1){ if(TimeStep == NbTimeStep - 1) { fprintf(PostStream, "};\n"); }
fprintf(PostStream, "};\n") ; break;
}
break ;
} }
} }
static void Gmsh_PrintElement(double Time, int TimeStep, int NbTimeStep, int NbH static void Gmsh_PrintElement(double Time, int TimeStep, int NbTimeStep,
armonic, int NbHarmonic, int HarmonicToTime, int Type,
int HarmonicToTime, int Type, int ElementNum, int int ElementNum, int NbrNodes, double *x,
NbrNodes, double *y, double *z, struct Value *Value,
double *x, double *y, double *z, struct Value *Val
ue,
struct PostSubOperation *PSO_P, int Store) struct PostSubOperation *PSO_P, int Store)
{ {
int i,j,k,jj ; int i, j, k, jj;
double TimeMH ; double TimeMH;
struct Value TmpValue ; struct Value TmpValue;
static std::vector<double> *Current_L ; static std::vector<double> *Current_L;
int symIndex[9] = {0, 1, 2, 1, 3, 4, 2, 4, 5} ; int symIndex[9] = {0, 1, 2, 1, 3, 4, 2, 4, 5};
int diagIndex[9] = {0, -1, -1, -1, 1, -1, -1, -1, 2} ; int diagIndex[9] = {0, -1, -1, -1, 1, -1, -1, -1, 2};
if(Gmsh_StartNewView){ if(Gmsh_StartNewView) {
Gmsh_StartNewView = 0 ; Gmsh_StartNewView = 0;
Gmsh_ResetStaticLists(); Gmsh_ResetStaticLists();
} }
switch (Value[0].Type) { switch(Value[0].Type) {
case SCALAR:
case SCALAR :
if(TimeStep == 0){ if(TimeStep == 0) {
switch(Type){ switch(Type) {
case POINT_ELEMENT : Current_L = &SP ; NbSP++ ; break ; case POINT_ELEMENT:
case LINE : Current_L = &SL ; NbSL++ ; break ; Current_L = &SP;
case TRIANGLE : Current_L = &ST ; NbST++ ; break ; NbSP++;
case QUADRANGLE : Current_L = &SQ ; NbSQ++ ; break ; break;
case TETRAHEDRON : Current_L = &SS ; NbSS++ ; break ; case LINE:
case HEXAHEDRON : Current_L = &SH ; NbSH++ ; break ; Current_L = &SL;
case PRISM : Current_L = &SI ; NbSI++ ; break ; NbSL++;
case PYRAMID : Current_L = &SY ; NbSY++ ; break ; break;
case TRIANGLE:
case LINE_2 : Current_L = &SL ; NbSL++ ; break ; Current_L = &ST;
case TRIANGLE_2 : Current_L = &ST ; NbST++ ; break ; NbST++;
case QUADRANGLE_2 : Current_L = &SQ ; NbSQ++ ; break ; break;
case QUADRANGLE_2_8N: Current_L = &SQ ; NbSQ++ ; break ; case QUADRANGLE:
case TETRAHEDRON_2 : Current_L = &SS ; NbSS++ ; break ; Current_L = &SQ;
case HEXAHEDRON_2 : Current_L = &SH ; NbSH++ ; break ; NbSQ++;
case HEXAHEDRON_2_20N: Current_L = &SH ; NbSH++ ; break ; break;
case PRISM_2 : Current_L = &SI ; NbSI++ ; break ; case TETRAHEDRON:
case PRISM_2_15N : Current_L = &SI ; NbSI++ ; break ; Current_L = &SS;
case PYRAMID_2 : Current_L = &SY ; NbSY++ ; break ; NbSS++;
case PYRAMID_2_13N : Current_L = &SY ; NbSY++ ; break ; break;
case HEXAHEDRON:
case LINE_3 : Current_L = &SL ; NbSL++ ; break ; Current_L = &SH;
case TRIANGLE_3 : Current_L = &ST ; NbST++ ; break ; NbSH++;
case QUADRANGLE_3 : Current_L = &SQ ; NbSQ++ ; break ; break;
case TETRAHEDRON_3 : Current_L = &SS ; NbSS++ ; break ; case PRISM:
case HEXAHEDRON_3 : Current_L = &SH ; NbSH++ ; break ; Current_L = &SI;
case PRISM_3 : Current_L = &SI ; NbSI++ ; break ; NbSI++;
case PYRAMID_3 : Current_L = &SY ; NbSY++ ; break ; break;
case PYRAMID:
case LINE_4 : Current_L = &SL ; NbSL++ ; break ; Current_L = &SY;
case TRIANGLE_4 : Current_L = &ST ; NbST++ ; break ; NbSY++;
case QUADRANGLE_4 : Current_L = &SQ ; NbSQ++ ; break ; break;
case TETRAHEDRON_4 : Current_L = &SS ; NbSS++ ; break ;
case HEXAHEDRON_4 : Current_L = &SH ; NbSH++ ; break ; case LINE_2:
case PRISM_4 : Current_L = &SI ; NbSI++ ; break ; Current_L = &SL;
//case PYRAMID_4 : Current_L = &SY ; NbSY++ ; break ; NbSL++;
} break;
if(Flag_GMSH_VERSION != 2){ case TRIANGLE_2:
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(x[i]); Current_L = &ST;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(y[i]); NbST++;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(z[i]); break;
case QUADRANGLE_2:
Current_L = &SQ;
NbSQ++;
break;
case QUADRANGLE_2_8N:
Current_L = &SQ;
NbSQ++;
break;
case TETRAHEDRON_2:
Current_L = &SS;
NbSS++;
break;
case HEXAHEDRON_2:
Current_L = &SH;
NbSH++;
break;
case HEXAHEDRON_2_20N:
Current_L = &SH;
NbSH++;
break;
case PRISM_2:
Current_L = &SI;
NbSI++;
break;
case PRISM_2_15N:
Current_L = &SI;
NbSI++;
break;
case PYRAMID_2:
Current_L = &SY;
NbSY++;
break;
case PYRAMID_2_13N:
Current_L = &SY;
NbSY++;
break;
case LINE_3:
Current_L = &SL;
NbSL++;
break;
case TRIANGLE_3:
Current_L = &ST;
NbST++;
break;
case QUADRANGLE_3:
Current_L = &SQ;
NbSQ++;
break;
case TETRAHEDRON_3:
Current_L = &SS;
NbSS++;
break;
case HEXAHEDRON_3:
Current_L = &SH;
NbSH++;
break;
case PRISM_3:
Current_L = &SI;
NbSI++;
break;
case PYRAMID_3:
Current_L = &SY;
NbSY++;
break;
case LINE_4:
Current_L = &SL;
NbSL++;
break;
case TRIANGLE_4:
Current_L = &ST;
NbST++;
break;
case QUADRANGLE_4:
Current_L = &SQ;
NbSQ++;
break;
case TETRAHEDRON_4:
Current_L = &SS;
NbSS++;
break;
case HEXAHEDRON_4:
Current_L = &SH;
NbSH++;
break;
case PRISM_4:
Current_L = &SI;
NbSI++;
break;
// case PYRAMID_4 : Current_L = &SY ; NbSY++ ; break ;
}
if(Flag_GMSH_VERSION != 2) {
for(i = 0; i < NbrNodes; i++) Current_L->push_back(x[i]);
for(i = 0; i < NbrNodes; i++) Current_L->push_back(y[i]);
for(i = 0; i < NbrNodes; i++) Current_L->push_back(z[i]);
} }
else{ else {
double tmp = ElementNum; double tmp = ElementNum;
Current_L->push_back(tmp); Current_L->push_back(tmp);
} }
} }
if (HarmonicToTime == 1) if(HarmonicToTime == 1)
for(k = 0 ; k < NbHarmonic ; k++){ for(k = 0; k < NbHarmonic; k++) {
List_Put(TimeValue_L, NbHarmonic*TimeStep+k, &Time); List_Put(TimeValue_L, NbHarmonic * TimeStep + k, &Time);
for(i = 0 ; i < NbrNodes ; i++) for(i = 0; i < NbrNodes; i++)
Current_L->push_back(Value[i].Val[MAX_DIM*k]); Current_L->push_back(Value[i].Val[MAX_DIM * k]);
} }
else else
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
List_Put(TimeValue_L, HarmonicToTime*TimeStep+k, &Time); List_Put(TimeValue_L, HarmonicToTime * TimeStep + k, &Time);
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
Current_L->push_back(TmpValue.Val[0]); Current_L->push_back(TmpValue.Val[0]);
} }
} }
break ; break;
case VECTOR : case VECTOR:
if(TimeStep == 0){ if(TimeStep == 0) {
switch(Type){ switch(Type) {
case POINT_ELEMENT : Current_L = &VP ; NbVP++ ; break ; case POINT_ELEMENT:
Current_L = &VP;
case LINE : Current_L = &VL ; NbVL++ ; break ; NbVP++;
case TRIANGLE : Current_L = &VT ; NbVT++ ; break ; break;
case QUADRANGLE : Current_L = &VQ ; NbVQ++ ; break ;
case TETRAHEDRON : Current_L = &VS ; NbVS++ ; break ; case LINE:
case HEXAHEDRON : Current_L = &VH ; NbVH++ ; break ; Current_L = &VL;
case PRISM : Current_L = &VI ; NbVI++ ; break ; NbVL++;
case PYRAMID : Current_L = &VY ; NbVY++ ; break ; break;
case TRIANGLE:
case LINE_2 : Current_L = &VL ; NbVL++ ; break ; Current_L = &VT;
case TRIANGLE_2 : Current_L = &VT ; NbVT++ ; break ; NbVT++;
case QUADRANGLE_2 : Current_L = &VQ ; NbVQ++ ; break ; break;
case QUADRANGLE_2_8N: Current_L = &VQ ; NbVQ++ ; break ; case QUADRANGLE:
case TETRAHEDRON_2 : Current_L = &VS ; NbVS++ ; break ; Current_L = &VQ;
case HEXAHEDRON_2 : Current_L = &VH ; NbVH++ ; break ; NbVQ++;
case HEXAHEDRON_2_20N: Current_L = &VH ; NbVH++ ; break ; break;
case PRISM_2 : Current_L = &VI ; NbVI++ ; break ; case TETRAHEDRON:
case PRISM_2_15N : Current_L = &VI ; NbVI++ ; break ; Current_L = &VS;
case PYRAMID_2 : Current_L = &VY ; NbVY++ ; break ; NbVS++;
case PYRAMID_2_13N : Current_L = &VY ; NbVY++ ; break ; break;
case HEXAHEDRON:
case LINE_3 : Current_L = &VL ; NbVL++ ; break ; Current_L = &VH;
case TRIANGLE_3 : Current_L = &VT ; NbVT++ ; break ; NbVH++;
case QUADRANGLE_3 : Current_L = &VQ ; NbVQ++ ; break ; break;
case TETRAHEDRON_3 : Current_L = &VS ; NbVS++ ; break ; case PRISM:
case HEXAHEDRON_3 : Current_L = &VH ; NbVH++ ; break ; Current_L = &VI;
case PRISM_3 : Current_L = &VI ; NbVI++ ; break ; NbVI++;
case PYRAMID_3 : Current_L = &VY ; NbVY++ ; break ; break;
case PYRAMID:
case LINE_4 : Current_L = &VL ; NbVL++ ; break ; Current_L = &VY;
case TRIANGLE_4 : Current_L = &VT ; NbVT++ ; break ; NbVY++;
case QUADRANGLE_4 : Current_L = &VQ ; NbVQ++ ; break ; break;
case TETRAHEDRON_4 : Current_L = &VS ; NbVS++ ; break ;
case HEXAHEDRON_4 : Current_L = &VH ; NbVH++ ; break ; case LINE_2:
case PRISM_4 : Current_L = &VI ; NbVI++ ; break ; Current_L = &VL;
//case PYRAMID_4 : Current_L = &VY ; NbVY++ ; break ; NbVL++;
} break;
if(Flag_GMSH_VERSION != 2){ case TRIANGLE_2:
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(x[i]); Current_L = &VT;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(y[i]); NbVT++;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(z[i]); break;
case QUADRANGLE_2:
Current_L = &VQ;
NbVQ++;
break;
case QUADRANGLE_2_8N:
Current_L = &VQ;
NbVQ++;
break;
case TETRAHEDRON_2:
Current_L = &VS;
NbVS++;
break;
case HEXAHEDRON_2:
Current_L = &VH;
NbVH++;
break;
case HEXAHEDRON_2_20N:
Current_L = &VH;
NbVH++;
break;
case PRISM_2:
Current_L = &VI;
NbVI++;
break;
case PRISM_2_15N:
Current_L = &VI;
NbVI++;
break;
case PYRAMID_2:
Current_L = &VY;
NbVY++;
break;
case PYRAMID_2_13N:
Current_L = &VY;
NbVY++;
break;
case LINE_3:
Current_L = &VL;
NbVL++;
break;
case TRIANGLE_3:
Current_L = &VT;
NbVT++;
break;
case QUADRANGLE_3:
Current_L = &VQ;
NbVQ++;
break;
case TETRAHEDRON_3:
Current_L = &VS;
NbVS++;
break;
case HEXAHEDRON_3:
Current_L = &VH;
NbVH++;
break;
case PRISM_3:
Current_L = &VI;
NbVI++;
break;
case PYRAMID_3:
Current_L = &VY;
NbVY++;
break;
case LINE_4:
Current_L = &VL;
NbVL++;
break;
case TRIANGLE_4:
Current_L = &VT;
NbVT++;
break;
case QUADRANGLE_4:
Current_L = &VQ;
NbVQ++;
break;
case TETRAHEDRON_4:
Current_L = &VS;
NbVS++;
break;
case HEXAHEDRON_4:
Current_L = &VH;
NbVH++;
break;
case PRISM_4:
Current_L = &VI;
NbVI++;
break;
// case PYRAMID_4 : Current_L = &VY ; NbVY++ ; break ;
}
if(Flag_GMSH_VERSION != 2) {
for(i = 0; i < NbrNodes; i++) Current_L->push_back(x[i]);
for(i = 0; i < NbrNodes; i++) Current_L->push_back(y[i]);
for(i = 0; i < NbrNodes; i++) Current_L->push_back(z[i]);
} }
else{ else {
double tmp = ElementNum; double tmp = ElementNum;
Current_L->push_back(tmp); Current_L->push_back(tmp);
} }
} }
if (HarmonicToTime == 1) if(HarmonicToTime == 1)
for(k = 0 ; k < NbHarmonic ; k++){ for(k = 0; k < NbHarmonic; k++) {
List_Put(TimeValue_L, NbHarmonic*TimeStep+k, &Time); List_Put(TimeValue_L, NbHarmonic * TimeStep + k, &Time);
for(i = 0 ; i < NbrNodes ; i++) for(i = 0; i < NbrNodes; i++)
for(j = 0 ; j < 3 ; j++) for(j = 0; j < 3; j++)
Current_L->push_back(Value[i].Val[MAX_DIM*k+j]); Current_L->push_back(Value[i].Val[MAX_DIM * k + j]);
} }
else else
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
List_Put(TimeValue_L, HarmonicToTime*TimeStep+k, &Time); List_Put(TimeValue_L, HarmonicToTime * TimeStep + k, &Time);
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
for(j = 0 ; j < 3 ; j++) for(j = 0; j < 3; j++) Current_L->push_back(TmpValue.Val[j]);
Current_L->push_back(TmpValue.Val[j]); }
} }
} break;
break ;
case TENSOR_DIAG:
case TENSOR_DIAG : case TENSOR_SYM:
case TENSOR_SYM : case TENSOR:
case TENSOR :
if(TimeStep == 0) {
if(TimeStep == 0){ switch(Type) {
switch(Type){ case POINT_ELEMENT:
case POINT_ELEMENT : Current_L = &TP ; NbTP++ ; break ; Current_L = &TP;
NbTP++;
case LINE : Current_L = &TL ; NbTL++ ; break ; break;
case TRIANGLE : Current_L = &TT ; NbTT++ ; break ;
case QUADRANGLE : Current_L = &TQ ; NbTQ++ ; break ; case LINE:
case TETRAHEDRON : Current_L = &TS1 ; NbTS++ ; break ; Current_L = &TL;
case HEXAHEDRON : Current_L = &TH ; NbTH++ ; break ; NbTL++;
case PRISM : Current_L = &TI ; NbTI++ ; break ; break;
case PYRAMID : Current_L = &TY ; NbTY++ ; break ; case TRIANGLE:
Current_L = &TT;
case LINE_2 : Current_L = &TL ; NbTL++ ; break ; NbTT++;
case TRIANGLE_2 : Current_L = &TT ; NbTT++ ; break ; break;
case QUADRANGLE_2 : Current_L = &TQ ; NbTQ++ ; break ; case QUADRANGLE:
case QUADRANGLE_2_8N: Current_L = &TQ ; NbTQ++ ; break ; Current_L = &TQ;
case TETRAHEDRON_2 : Current_L = &TS1 ; NbTS++ ; break ; NbTQ++;
case HEXAHEDRON_2 : Current_L = &TH ; NbTH++ ; break ; break;
case HEXAHEDRON_2_20N: Current_L = &TH ; NbTH++ ; break ; case TETRAHEDRON:
case PRISM_2 : Current_L = &TI ; NbTI++ ; break ; Current_L = &TS1;
case PRISM_2_15N : Current_L = &TI ; NbTI++ ; break ; NbTS++;
case PYRAMID_2 : Current_L = &TY ; NbTY++ ; break ; break;
case PYRAMID_2_13N : Current_L = &TY ; NbTY++ ; break ; case HEXAHEDRON:
Current_L = &TH;
case LINE_3 : Current_L = &TL ; NbTL++ ; break ; NbTH++;
case TRIANGLE_3 : Current_L = &TT ; NbTT++ ; break ; break;
case QUADRANGLE_3 : Current_L = &TQ ; NbTQ++ ; break ; case PRISM:
case TETRAHEDRON_3 : Current_L = &TS1 ; NbTS++ ; break ; Current_L = &TI;
case HEXAHEDRON_3 : Current_L = &TH ; NbTH++ ; break ; NbTI++;
case PRISM_3 : Current_L = &TI ; NbTI++ ; break ; break;
case PYRAMID_3 : Current_L = &TY ; NbTY++ ; break ; case PYRAMID:
Current_L = &TY;
case LINE_4 : Current_L = &TL ; NbTL++ ; break ; NbTY++;
case TRIANGLE_4 : Current_L = &TT ; NbTT++ ; break ; break;
case QUADRANGLE_4 : Current_L = &TQ ; NbTQ++ ; break ;
case TETRAHEDRON_4 : Current_L = &TS1 ; NbTS++ ; break ; case LINE_2:
case HEXAHEDRON_4 : Current_L = &TH ; NbTH++ ; break ; Current_L = &TL;
case PRISM_4 : Current_L = &TI ; NbTI++ ; break ; NbTL++;
//case PYRAMID_4 : Current_L = &TY ; NbTY++ ; break ; break;
} case TRIANGLE_2:
if(Flag_GMSH_VERSION != 2){ Current_L = &TT;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(x[i]); NbTT++;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(y[i]); break;
for(i = 0 ; i < NbrNodes ; i++) Current_L->push_back(z[i]); case QUADRANGLE_2:
Current_L = &TQ;
NbTQ++;
break;
case QUADRANGLE_2_8N:
Current_L = &TQ;
NbTQ++;
break;
case TETRAHEDRON_2:
Current_L = &TS1;
NbTS++;
break;
case HEXAHEDRON_2:
Current_L = &TH;
NbTH++;
break;
case HEXAHEDRON_2_20N:
Current_L = &TH;
NbTH++;
break;
case PRISM_2:
Current_L = &TI;
NbTI++;
break;
case PRISM_2_15N:
Current_L = &TI;
NbTI++;
break;
case PYRAMID_2:
Current_L = &TY;
NbTY++;
break;
case PYRAMID_2_13N:
Current_L = &TY;
NbTY++;
break;
case LINE_3:
Current_L = &TL;
NbTL++;
break;
case TRIANGLE_3:
Current_L = &TT;
NbTT++;
break;
case QUADRANGLE_3:
Current_L = &TQ;
NbTQ++;
break;
case TETRAHEDRON_3:
Current_L = &TS1;
NbTS++;
break;
case HEXAHEDRON_3:
Current_L = &TH;
NbTH++;
break;
case PRISM_3:
Current_L = &TI;
NbTI++;
break;
case PYRAMID_3:
Current_L = &TY;
NbTY++;
break;
case LINE_4:
Current_L = &TL;
NbTL++;
break;
case TRIANGLE_4:
Current_L = &TT;
NbTT++;
break;
case QUADRANGLE_4:
Current_L = &TQ;
NbTQ++;
break;
case TETRAHEDRON_4:
Current_L = &TS1;
NbTS++;
break;
case HEXAHEDRON_4:
Current_L = &TH;
NbTH++;
break;
case PRISM_4:
Current_L = &TI;
NbTI++;
break;
// case PYRAMID_4 : Current_L = &TY ; NbTY++ ; break ;
}
if(Flag_GMSH_VERSION != 2) {
for(i = 0; i < NbrNodes; i++) Current_L->push_back(x[i]);
for(i = 0; i < NbrNodes; i++) Current_L->push_back(y[i]);
for(i = 0; i < NbrNodes; i++) Current_L->push_back(z[i]);
} }
else{ else {
double tmp = ElementNum; double tmp = ElementNum;
Current_L->push_back(tmp); Current_L->push_back(tmp);
} }
} }
if (HarmonicToTime == 1) if(HarmonicToTime == 1)
for(k = 0 ; k < NbHarmonic ; k++){ for(k = 0; k < NbHarmonic; k++) {
List_Put(TimeValue_L, NbHarmonic*TimeStep+k, &Time); List_Put(TimeValue_L, NbHarmonic * TimeStep + k, &Time);
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
for(j = 0 ; j < 9 ; j++){ for(j = 0; j < 9; j++) {
if(Value[0].Type != TENSOR_DIAG) { if(Value[0].Type != TENSOR_DIAG) {
if(Value[0].Type == TENSOR_SYM) jj = symIndex[j]; if(Value[0].Type == TENSOR_SYM)
else jj = j; jj = symIndex[j];
Current_L->push_back(Value[i].Val[MAX_DIM*k+jj]); else
jj = j;
Current_L->push_back(Value[i].Val[MAX_DIM * k + jj]);
} }
else { else {
jj = diagIndex[j]; jj = diagIndex[j];
if(jj == -1) Current_L->push_back(0.); if(jj == -1)
else Current_L->push_back(Value[i].Val[MAX_DIM*k+jj]); Current_L->push_back(0.);
else
Current_L->push_back(Value[i].Val[MAX_DIM * k + jj]);
} }
} }
} }
} }
else else
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
List_Put(TimeValue_L, HarmonicToTime*TimeStep+k, &Time); List_Put(TimeValue_L, HarmonicToTime * TimeStep + k, &Time);
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
for(j = 0 ; j < 9 ; j++) { for(j = 0; j < 9; j++) {
if(Value[0].Type != TENSOR_DIAG) { if(Value[0].Type != TENSOR_DIAG) {
if(Value[0].Type == TENSOR_SYM) jj = symIndex[j]; if(Value[0].Type == TENSOR_SYM)
else jj = j; jj = symIndex[j];
else
jj = j;
Current_L->push_back(TmpValue.Val[jj]); Current_L->push_back(TmpValue.Val[jj]);
} }
else { else {
jj = diagIndex[j]; jj = diagIndex[j];
if(jj == -1) Current_L->push_back(0.); if(jj == -1)
else Current_L->push_back(TmpValue.Val[jj]); Current_L->push_back(0.);
else
Current_L->push_back(TmpValue.Val[jj]);
} }
} }
} }
} }
} }
// reduce memory requirements by automatically partitioning large // reduce memory requirements by automatically partitioning large
// output views into chunks not larger than 1Gb // output views into chunks not larger than 1Gb
if(Flag_GMSH_VERSION == 2 && TimeStep == NbTimeStep - 1 && if(Flag_GMSH_VERSION == 2 && TimeStep == NbTimeStep - 1 &&
Current_L->size() > (int)(1024 * 1024 * 1024 / sizeof(double))){ Current_L->size() > (int)(1024 * 1024 * 1024 / sizeof(double))) {
Format_PostFooter(PSO_P, Store); Format_PostFooter(PSO_P, Store);
CurrentPartitionNumber++; CurrentPartitionNumber++;
Gmsh_StartNewView = 1; Gmsh_StartNewView = 1;
} }
} }
static void dVecWrite(std::vector<double> &v, FILE *fp, bool binary) static void dVecWrite(std::vector<double> &v, FILE *fp, bool binary)
{ {
if(v.empty()) return; if(v.empty()) return;
if(binary) if(binary)
fwrite(&v[0], sizeof(double), v.size(), fp); fwrite(&v[0], sizeof(double), v.size(), fp);
else else
for(unsigned i = 0; i < v.size(); i++) for(unsigned i = 0; i < v.size(); i++) fprintf(fp, " %.16g", v[i]);
fprintf(fp, " %.16g", v[i]);
} }
static void cVecWrite(std::vector<char> &v, FILE *fp, bool binary) static void cVecWrite(std::vector<char> &v, FILE *fp, bool binary)
{ {
if(v.empty()) return; if(v.empty()) return;
if(binary) if(binary)
fwrite(&v[0], sizeof(char), v.size(), fp); fwrite(&v[0], sizeof(char), v.size(), fp);
else else
for(unsigned i = 0; i < v.size(); i++) for(unsigned i = 0; i < v.size(); i++) fputc(v[i], fp);
fputc(v[i], fp);
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Gnuplot format */ /* Gnuplot format */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
static void Gnuplot_PrintElement(int Format, double Time, int TimeStep, int NbrT static void Gnuplot_PrintElement(int Format, double Time, int TimeStep,
imeSteps, int NbrTimeSteps, int NbrHarmonics,
int NbrHarmonics, int HarmonicToTime, int HarmonicToTime, int ElementType,
int ElementType, int NumElement, int NbrNodes, int NumElement, int NbrNodes, double *x,
double *x, double *y, double *z, double *Dummy, double *y, double *z, double *Dummy,
struct Value *Value) struct Value *Value)
{ {
static int Size, TmpIndex ; static int Size, TmpIndex;
static double * TmpValues ; static double *TmpValues;
int i, j, k, t, i2, k2 ; int i, j, k, t, i2, k2;
double TimeMH ; double TimeMH;
struct Value TmpValue ; struct Value TmpValue;
if(!PostStream) return; if(!PostStream) return;
if(TimeStep == 0){ if(TimeStep == 0) {
switch(Value->Type){ switch(Value->Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
TmpValues = (double*) Malloc(NbrTimeSteps*NbrNodes*NbrHarmonics*Size*sizeof( TmpValues = (double *)Malloc(NbrTimeSteps * NbrNodes * NbrHarmonics * Size *
double)); sizeof(double));
TmpIndex = 0; TmpIndex = 0;
} }
for(i = 0 ; i < NbrNodes ; i++) for(i = 0; i < NbrNodes; i++)
for(k = 0 ; k < NbrHarmonics ; k++) for(k = 0; k < NbrHarmonics; k++)
for(j = 0 ; j < Size ; j++) for(j = 0; j < Size; j++)
TmpValues[TmpIndex++] = Value[i].Val[MAX_DIM*k+j]; TmpValues[TmpIndex++] = Value[i].Val[MAX_DIM * k + j];
if(TimeStep == NbrTimeSteps-1){ if(TimeStep == NbrTimeSteps - 1) {
for(i = 0; i <= NbrNodes;
for(i = 0 ; i <= NbrNodes ; i++){ /* New line for each node, closed loop for i++) { /* New line for each node, closed loop for tri/qua */
tri/qua */
if(i != NbrNodes) if(i != NbrNodes)
i2 = i ; i2 = i;
else{ else {
if(NbrNodes < 3) break ; if(NbrNodes < 3)
else i2 = 0 ; break;
else
i2 = 0;
} }
fprintf(PostStream, "%d %d ", GetDP2Gmsh(ElementType), NumElement); fprintf(PostStream, "%d %d ", GetDP2Gmsh(ElementType), NumElement);
fprintf(PostStream, " %.16g %.16g %.16g ", x[i2], y[i2], z[i2]); fprintf(PostStream, " %.16g %.16g %.16g ", x[i2], y[i2], z[i2]);
if(Dummy){ if(Dummy) {
if(Dummy[3]<0){ if(Dummy[3] < 0) {
if(!i) if(!i)
fprintf(PostStream, " %.16g %.16g 0 ", Dummy[0], Dummy[2]); fprintf(PostStream, " %.16g %.16g 0 ", Dummy[0], Dummy[2]);
else else
fprintf(PostStream, " %.16g %.16g 0 ", Dummy[1], Dummy[2]); fprintf(PostStream, " %.16g %.16g 0 ", Dummy[1], Dummy[2]);
} }
else else
fprintf(PostStream, " %.16g %.16g %.16g ", Dummy[0], Dummy[1], Dummy[2] fprintf(PostStream, " %.16g %.16g %.16g ", Dummy[0], Dummy[1],
); Dummy[2]);
} }
else else
fprintf(PostStream, " 0 0 0 "); fprintf(PostStream, " 0 0 0 ");
for(t = 0 ; t < NbrTimeSteps ; t++){ for(t = 0; t < NbrTimeSteps; t++) {
if(HarmonicToTime == 1) {
if (HarmonicToTime == 1) { for(k = 0; k < NbrHarmonics; k++) {
for(k = 0 ; k < NbrHarmonics ; k++) { for(j = 0; j < Size; j++) {
for(j = 0 ; j < Size ; j++){ fprintf(PostStream, " %.16g",
fprintf(PostStream, " %.16g", TmpValues[t * NbrNodes * NbrHarmonics * Size +
TmpValues[ t*NbrNodes*NbrHarmonics*Size i2 * NbrHarmonics * Size + k * Size + j]);
+ i2*NbrHarmonics*Size }
+ k*Size fprintf(PostStream, " ");
+ j ]); }
} }
fprintf(PostStream, " "); else {
} TmpValue.Type = Value->Type;
} for(k = 0; k < HarmonicToTime; k++) {
else { for(k2 = 0; k2 < NbrHarmonics; k2++)
TmpValue.Type = Value->Type ; for(j = 0; j < Size; j++)
for(k = 0 ; k < HarmonicToTime ; k++){ TmpValue.Val[MAX_DIM * k2 + j] =
TmpValues[t * NbrNodes * NbrHarmonics * Size +
for(k2 = 0 ; k2 < NbrHarmonics ; k2++) i2 * NbrHarmonics * Size + k2 * Size + j];
for(j = 0 ; j < Size ; j++)
TmpValue.Val[MAX_DIM*k2+j] = F_MHToTime0(k, &TmpValue, &TmpValue, k, HarmonicToTime, &TimeMH);
TmpValues[ t*NbrNodes*NbrHarmonics*Size for(j = 0; j < Size; j++)
+ i2*NbrHarmonics*Size fprintf(PostStream, "%.16g", TmpValue.Val[0]);
+ k2*Size fprintf(PostStream, " ");
+ j ] ; }
}
F_MHToTime0(k, &TmpValue, &TmpValue, k, HarmonicToTime, &TimeMH) ; fprintf(PostStream, " ");
for(j = 0 ; j < Size ; j++)
fprintf(PostStream, "%.16g", TmpValue.Val[0]);
fprintf(PostStream, " ");
}
}
fprintf(PostStream, " ");
} /* for t */ } /* for t */
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
} /* for i */ } /* for i */
if(NbrNodes > 1) fprintf(PostStream, "\n"); if(NbrNodes > 1) fprintf(PostStream, "\n");
Free(TmpValues); Free(TmpValues);
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* Tabular format */ /* Tabular format */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
// global static list for tables output // global static list for tables output
static int TableList_StartNew = 0; static int TableList_StartNew = 0;
static std::list<double> TableList; static std::list<double> TableList;
static void Tabular_PrintElement(struct PostSubOperation *PSO_P, static void Tabular_PrintElement(struct PostSubOperation *PSO_P, int Format,
int Format, double Time, int TimeStep, int NbrT double Time, int TimeStep, int NbrTimeSteps,
imeSteps,
int NbrHarmonics, int HarmonicToTime, int NbrHarmonics, int HarmonicToTime,
int ElementType, int NumElement, int NbrNodes, int ElementType, int NumElement, int NbrNodes,
double *x, double *y, double *z, double *Dummy, double *x, double *y, double *z, double *Dummy,
struct Value *Value) struct Value *Value)
{ {
static int Size ; static int Size;
int i,j,k ; int i, j, k;
double TimeMH ; double TimeMH;
struct Value TmpValue ; struct Value TmpValue;
if(TableList_StartNew){ if(TableList_StartNew) {
TableList_StartNew = 0 ; TableList_StartNew = 0;
TableList.clear(); TableList.clear();
} }
if(!PostStream) return; if(!PostStream) return;
if(TimeStep == 0){ if(TimeStep == 0) {
switch(Value->Type){ switch(Value->Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
} }
if(Format == FORMAT_SPACE_TABLE || Format == FORMAT_SIMPLE_SPACE_TABLE if(Format == FORMAT_SPACE_TABLE || Format == FORMAT_SIMPLE_SPACE_TABLE ||
|| Format == FORMAT_VALUE_ONLY){ Format == FORMAT_VALUE_ONLY) {
if(TimeStep == 0){ if(TimeStep == 0) {
if(Format != FORMAT_SIMPLE_SPACE_TABLE && Format != FORMAT_VALUE_ONLY){ if(Format != FORMAT_SIMPLE_SPACE_TABLE && Format != FORMAT_VALUE_ONLY) {
fprintf(PostStream, "%d %d ", GetDP2Gmsh(ElementType), NumElement); fprintf(PostStream, "%d %d ", GetDP2Gmsh(ElementType), NumElement);
TableList.push_back(GetDP2Gmsh(ElementType)); TableList.push_back(GetDP2Gmsh(ElementType));
TableList.push_back(NumElement); TableList.push_back(NumElement);
} }
if(Format != FORMAT_VALUE_ONLY) if(Format != FORMAT_VALUE_ONLY)
for(i=0 ; i<NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
fprintf(PostStream, "%.16g %.16g %.16g ", x[i], y[i], z[i]); fprintf(PostStream, "%.16g %.16g %.16g ", x[i], y[i], z[i]);
TableList.push_back(x[i]); TableList.push_back(x[i]);
TableList.push_back(y[i]); TableList.push_back(y[i]);
TableList.push_back(z[i]); TableList.push_back(z[i]);
} }
if(Format != FORMAT_SIMPLE_SPACE_TABLE && Format != FORMAT_VALUE_ONLY){ if(Format != FORMAT_SIMPLE_SPACE_TABLE && Format != FORMAT_VALUE_ONLY) {
if(Dummy){ if(Dummy) {
fprintf(PostStream, "%.16g %.16g %.16g ", Dummy[0], Dummy[1], Dummy[ fprintf(PostStream, "%.16g %.16g %.16g ", Dummy[0], Dummy[1],
2]); Dummy[2]);
TableList.push_back(Dummy[0]); TableList.push_back(Dummy[0]);
TableList.push_back(Dummy[1]); TableList.push_back(Dummy[1]);
TableList.push_back(Dummy[2]); TableList.push_back(Dummy[2]);
} }
else{ else {
fprintf(PostStream, "0 0 0 "); fprintf(PostStream, "0 0 0 ");
TableList.push_back(0); TableList.push_back(0);
TableList.push_back(0); TableList.push_back(0);
TableList.push_back(0); TableList.push_back(0);
} }
} }
} }
} }
if (HarmonicToTime == 1) { if(HarmonicToTime == 1) {
if(Format == FORMAT_TIME_TABLE){ if(Format == FORMAT_TIME_TABLE) {
fprintf(PostStream, "%d %.16g ", TimeStep, Time); fprintf(PostStream, "%d %.16g ", TimeStep, Time);
TableList.push_back(TimeStep); TableList.push_back(TimeStep);
TableList.push_back(Time); TableList.push_back(Time);
for(i=0 ; i<NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
fprintf(PostStream, " %.16g %.16g %.16g ", x[i], y[i], z[i]); fprintf(PostStream, " %.16g %.16g %.16g ", x[i], y[i], z[i]);
TableList.push_back(x[i]); TableList.push_back(x[i]);
TableList.push_back(y[i]); TableList.push_back(y[i]);
TableList.push_back(z[i]); TableList.push_back(z[i]);
} }
} }
for(k = 0 ; k < NbrHarmonics ; k++) { for(k = 0; k < NbrHarmonics; k++) {
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
for(j = 0 ; j < Size ; j++){ for(j = 0; j < Size; j++) {
if (Format != FORMAT_VALUE_ONLY){ if(Format != FORMAT_VALUE_ONLY) {
fprintf(PostStream, " %.16g", Value[i].Val[MAX_DIM*k+j]); fprintf(PostStream, " %.16g", Value[i].Val[MAX_DIM * k + j]);
TableList.push_back(Value[i].Val[MAX_DIM*k+j]); TableList.push_back(Value[i].Val[MAX_DIM * k + j]);
} }
else{ else {
fprintf(PostStream, " %s_%d_%d = %.16g", fprintf(PostStream, " %s_%d_%d = %.16g", PSO_P->ValueName, j,
PSO_P->ValueName, j, PSO_P->ValueIndex, Value[i].Val[MAX_DIM PSO_P->ValueIndex, Value[i].Val[MAX_DIM * k + j]);
*k+j]); TableList.push_back(Value[i].Val[MAX_DIM * k + j]);
TableList.push_back(Value[i].Val[MAX_DIM*k+j]); if(j < Size - 1) fprintf(PostStream, "\n");
if (j<Size-1)
fprintf(PostStream, "\n");
} }
} }
fprintf(PostStream, " "); fprintf(PostStream, " ");
} }
fprintf(PostStream, " "); fprintf(PostStream, " ");
} }
} }
else { else {
for(k = 0 ; k < HarmonicToTime ; k++){ for(k = 0; k < HarmonicToTime; k++) {
for(i = 0 ; i < NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
F_MHToTime0(k+i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH) ; F_MHToTime0(k + i, &Value[i], &TmpValue, k, HarmonicToTime, &TimeMH);
if (!i && Format == FORMAT_TIME_TABLE) { if(!i && Format == FORMAT_TIME_TABLE) {
fprintf(PostStream, "%d %.16g ", TimeStep, TimeMH); fprintf(PostStream, "%d %.16g ", TimeStep, TimeMH);
TableList.push_back(TimeStep); TableList.push_back(TimeStep);
TableList.push_back(TimeMH); TableList.push_back(TimeMH);
for(i=0 ; i<NbrNodes ; i++){ for(i = 0; i < NbrNodes; i++) {
fprintf(PostStream, " %.16g %.16g %.16g ", x[i], y[i], z[i]); fprintf(PostStream, " %.16g %.16g %.16g ", x[i], y[i], z[i]);
TableList.push_back(x[i]); TableList.push_back(x[i]);
TableList.push_back(y[i]); TableList.push_back(y[i]);
TableList.push_back(z[i]); TableList.push_back(z[i]);
} }
} }
for(j = 0 ; j < Size ; j++){ for(j = 0; j < Size; j++) {
fprintf(PostStream, " %.16g", TmpValue.Val[j]) ; fprintf(PostStream, " %.16g", TmpValue.Val[j]);
TableList.push_back(TmpValue.Val[j]); TableList.push_back(TmpValue.Val[j]);
} }
fprintf(PostStream, " "); fprintf(PostStream, " ");
} }
if(Format == FORMAT_TIME_TABLE) if(Format == FORMAT_TIME_TABLE) fprintf(PostStream, "\n");
fprintf(PostStream, "\n");
} }
} }
if(TimeStep == NbrTimeSteps-1 || Format == FORMAT_TIME_TABLE) if(TimeStep == NbrTimeSteps - 1 || Format == FORMAT_TIME_TABLE)
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* NodeTable format */ /* NodeTable format */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
// global static map for node table output (cannot be saved incrementally for // global static map for node table output (cannot be saved incrementally for
// each element) // each element)
static int NodeTable_StartNew = 0; static int NodeTable_StartNew = 0;
static std::map<int, std::vector<double> > NodeTable; static std::map<int, std::vector<double> > NodeTable;
static void NodeTable_PrintElement(int TimeStep, int NbTimeStep, int NbrHarmonic static void NodeTable_PrintElement(int TimeStep, int NbTimeStep,
s, int NbrHarmonics, struct PostElement *PE)
struct PostElement *PE)
{ {
if(NodeTable_StartNew){ if(NodeTable_StartNew) {
NodeTable_StartNew = 0 ; NodeTable_StartNew = 0;
NodeTable.clear(); NodeTable.clear();
} }
for(int i = 0 ; i < PE->NbrNodes ; i++){ for(int i = 0; i < PE->NbrNodes; i++) {
int n = PE->NumNodes[i]; int n = PE->NumNodes[i];
int Size = 0; int Size = 0;
switch(PE->Value[0].Type){ switch(PE->Value[0].Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
if(n > 0 && Size){ // we have data on an actual node if(n > 0 && Size) { // we have data on an actual node
NodeTable[n].resize(NbTimeStep * NbrHarmonics * Size, 0.); NodeTable[n].resize(NbTimeStep * NbrHarmonics * Size, 0.);
for(int k = 0 ; k < NbrHarmonics ; k++){ for(int k = 0; k < NbrHarmonics; k++) {
for(int j = 0 ; j < Size ; j++){ for(int j = 0; j < Size; j++) {
double val = PE->Value[i].Val[MAX_DIM * k + j]; double val = PE->Value[i].Val[MAX_DIM * k + j];
int idx = NbrHarmonics * Size * TimeStep + k * Size + j; int idx = NbrHarmonics * Size * TimeStep + k * Size + j;
NodeTable[n][idx] = val; NodeTable[n][idx] = val;
} }
} }
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* ElementTable format * / /* ElementTable format */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
// global static map for element table output (cannot be saved incrementally for // global static map for element table output (cannot be saved incrementally for
// each element) // each element)
static int ElementTable_StartNew = 0; static int ElementTable_StartNew = 0;
static std::map<int, std::vector<double> > ElementTable; static std::map<int, std::vector<double> > ElementTable;
static void ElementTable_PrintElement(int TimeStep, int NbTimeStep, int NbrHarmo static void ElementTable_PrintElement(int TimeStep, int NbTimeStep,
nics, int NbrHarmonics, struct PostElement *PE)
struct PostElement *PE)
{ {
if(ElementTable_StartNew){ if(ElementTable_StartNew) {
ElementTable_StartNew = 0 ; ElementTable_StartNew = 0;
ElementTable.clear(); ElementTable.clear();
} }
int numEle = -1; int numEle = -1;
if(PE->Index >= 0 && PE->Index < Geo_GetNbrGeoElements()){ if(PE->Index >= 0 && PE->Index < Geo_GetNbrGeoElements()) {
numEle = Geo_GetGeoElement(PE->Index)->Num; numEle = Geo_GetGeoElement(PE->Index)->Num;
int Size = 0; int Size = 0;
switch(PE->Value[0].Type){ switch(PE->Value[0].Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
if(Size){ if(Size) {
ElementTable[numEle].resize ElementTable[numEle].resize(
(NbTimeStep * PE->NbrNodes * NbrHarmonics * Size, 0.); NbTimeStep * PE->NbrNodes * NbrHarmonics * Size, 0.);
for(int i = 0 ; i < PE->NbrNodes ; i++){ for(int i = 0; i < PE->NbrNodes; i++) {
for(int k = 0 ; k < NbrHarmonics ; k++){ for(int k = 0; k < NbrHarmonics; k++) {
for(int j = 0 ; j < Size ; j++){ for(int j = 0; j < Size; j++) {
double val = PE->Value[i].Val[MAX_DIM * k + j]; double val = PE->Value[i].Val[MAX_DIM * k + j];
int idx = TimeStep * PE->NbrNodes * NbrHarmonics * Size + int idx = TimeStep * PE->NbrNodes * NbrHarmonics * Size +
i * NbrHarmonics * Size + k * Size + j; i * NbrHarmonics * Size + k * Size + j;
ElementTable[numEle][idx] = val; ElementTable[numEle][idx] = val;
} }
} }
} }
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* S t o r e P o s t O p R e s u l t */ /* S t o r e P o s t O p R e s u l t */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
static List_T *PostOpResults_L=NULL; static List_T *PostOpResults_L = NULL;
static void StorePostOpResult(int NbrHarmonics, struct PostElement *PE) static void StorePostOpResult(int NbrHarmonics, struct PostElement *PE)
{ {
int Size; int Size;
double val; double val;
if(!PostOpResults_L) if(!PostOpResults_L)
PostOpResults_L = List_Create(1000,1000,sizeof(double)); PostOpResults_L = List_Create(1000, 1000, sizeof(double));
for(int i = 0 ; i < PE->NbrNodes ; i++){ for(int i = 0; i < PE->NbrNodes; i++) {
Size = 0; Size = 0;
switch(PE->Value[0].Type){ switch(PE->Value[0].Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
if(Size){ // we have data if(Size) { // we have data
for(int k = 0 ; k < NbrHarmonics ; k++){ for(int k = 0; k < NbrHarmonics; k++) {
for(int j = 0 ; j < Size ; j++){ for(int j = 0; j < Size; j++) {
val = PE->Value[i].Val[MAX_DIM * k + j]; val = PE->Value[i].Val[MAX_DIM * k + j];
List_Add(PostOpResults_L, &val); List_Add(PostOpResults_L, &val);
} }
} }
} }
} }
} }
static void StorePostOpResult(int NbrHarmonics, struct Value *Value) static void StorePostOpResult(int NbrHarmonics, struct Value *Value)
{ {
int Size; int Size;
double val; double val;
if(!PostOpResults_L) if(!PostOpResults_L)
PostOpResults_L = List_Create(1000,1000,sizeof(double)); PostOpResults_L = List_Create(1000, 1000, sizeof(double));
Size = 0; Size = 0;
switch(Value[0].Type){ switch(Value[0].Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
if(Size){ // we have data if(Size) { // we have data
for(int k = 0 ; k < NbrHarmonics ; k++){ for(int k = 0; k < NbrHarmonics; k++) {
for(int j = 0 ; j < Size ; j++){ for(int j = 0; j < Size; j++) {
val = Value[0].Val[MAX_DIM * k + j]; val = Value[0].Val[MAX_DIM * k + j];
List_Add(PostOpResults_L, &val); List_Add(PostOpResults_L, &val);
} }
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* UNV format */ /* UNV format */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
#if !defined(HAVE_NX) #if !defined(HAVE_NX)
#define NX { Message::Error("UNV export not available in this version"); } #define NX \
{ \
Message::Error("UNV export not available in this version"); \
}
#else #else
#define NX ; #define NX ;
#endif #endif
double NXUnv_UnitFactor = 1; double NXUnv_UnitFactor = 1;
int NXUnv_DatasetLocation = 3; //1: Data at nodes, 2: Data on elements, 3: Da int NXUnv_DatasetLocation =
ta at nodes on elements 3; // 1: Data at nodes, 2: Data on elements, 3: Data at nodes on elements
void Unv_PrintHeader(FILE *PostStream, char *name, double Time, int TimeStep, do void Unv_PrintHeader(FILE *PostStream, char *name, double Time, int TimeStep,
uble& NXUnv_UnitFactor, int& NXUnv_DatasetLocation) NX double &NXUnv_UnitFactor, int &NXUnv_DatasetLocation) NX
void Unv_PrintFooter(FILE *PostStream) NX void Unv_PrintFooter(FILE *PostStream) NX
void Unv_PrintElement(FILE *PostStream, int Num_Element, int NbrNodes, struct Va void Unv_PrintElement(FILE *PostStream, int Num_Element, int NbrNodes,
lue *Value, int NbrHarmonics, int& NXUnv_DatasetLocation, double& NXUnv_UnitFact struct Value *Value, int NbrHarmonics,
or) NX int &NXUnv_DatasetLocation, double &NXUnv_UnitFactor) NX
void Unv_PrintNodeTable(FILE *PostStream, std::map<int, std::vector<double>> &No void Unv_PrintNodeTable(FILE *PostStream,
deTable, double& NXUnv_UnitFactor) NX std::map<int, std::vector<double> > &NodeTable,
void Unv_PrintRegion(FILE *PostStream, int Flag_Comma, int numRegion, int NbrHar double &NXUnv_UnitFactor) NX
monics, int Size, struct Value *Value, double& NXUnv_UnitFactor) NX void Unv_PrintRegion(FILE *PostStream, int Flag_Comma, int numRegion,
int NbrHarmonics, int Size, struct Value *Value,
double &NXUnv_UnitFactor) NX
#undef NX #undef NX
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* F o r m a t _ P o s t F o r m a t */ /* F o r m a t _ P o s t F o r m a t */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Format_PostFormat(struct PostSubOperation *PSO_P) void Format_PostFormat(struct PostSubOperation *PSO_P)
{ {
if(!PostStream || PSO_P->Type == POP_EXPRESSION) return; if(!PostStream || PSO_P->Type == POP_EXPRESSION) return;
int Format = PSO_P->Format; int Format = PSO_P->Format;
int NoMesh = PSO_P->NoMesh; int NoMesh = PSO_P->NoMesh;
switch(Format){ switch(Format) {
case FORMAT_GMSH : case FORMAT_GMSH:
if((PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0) && if((PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0) &&
!PSO_P->FileOut) break; !PSO_P->FileOut)
if(Flag_GMSH_VERSION == 2){ break;
fprintf(PostStream, "$MeshFormat\n") ; if(Flag_GMSH_VERSION == 2) {
fprintf(PostStream, "2.2 %d %d\n", Flag_BIN, (int)sizeof(double)) ; fprintf(PostStream, "$MeshFormat\n");
if(Flag_BIN){ fprintf(PostStream, "2.2 %d %d\n", Flag_BIN, (int)sizeof(double));
int one=1; if(Flag_BIN) {
int one = 1;
fwrite(&one, sizeof(int), 1, PostStream); fwrite(&one, sizeof(int), 1, PostStream);
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
} }
fprintf(PostStream, "$EndMeshFormat\n") ; fprintf(PostStream, "$EndMeshFormat\n");
if(!NoMesh){ if(!NoMesh) {
std::vector<Geo_Element*> elements; std::vector<Geo_Element *> elements;
std::set<int> nodes; std::set<int> nodes;
if(PSO_P->SubType == PRINT_ONELEMENTSOF){ if(PSO_P->SubType == PRINT_ONELEMENTSOF) {
List_T *Region_L = List_T *Region_L =
((struct Group *)List_Pointer(Problem_S.Group, ((struct Group *)List_Pointer(Problem_S.Group,
PSO_P->Case.OnRegion.RegionIndex))->In PSO_P->Case.OnRegion.RegionIndex))
itialList ; ->InitialList;
for(int i = 0 ; i < Geo_GetNbrGeoElements() ; i++) { for(int i = 0; i < Geo_GetNbrGeoElements(); i++) {
Geo_Element *Geo_Element = Geo_GetGeoElement(i) ; Geo_Element *Geo_Element = Geo_GetGeoElement(i);
if(List_Search(Region_L, &Geo_Element->Region, fcmp_int)){ if(List_Search(Region_L, &Geo_Element->Region, fcmp_int)) {
elements.push_back(Geo_Element); elements.push_back(Geo_Element);
for (int j = 0 ; j < Geo_Element->NbrNodes ; j++) for(int j = 0; j < Geo_Element->NbrNodes; j++)
nodes.insert(Geo_Element->NumNodes[j]) ; nodes.insert(Geo_Element->NumNodes[j]);
} }
} }
} }
else{ else {
for(int i = 0 ; i < Geo_GetNbrGeoElements() ; i++) { for(int i = 0; i < Geo_GetNbrGeoElements(); i++) {
Geo_Element *Geo_Element = Geo_GetGeoElement(i) ; Geo_Element *Geo_Element = Geo_GetGeoElement(i);
elements.push_back(Geo_Element); elements.push_back(Geo_Element);
for (int j = 0 ; j < Geo_Element->NbrNodes ; j++) for(int j = 0; j < Geo_Element->NbrNodes; j++)
nodes.insert(Geo_Element->NumNodes[j]) ; nodes.insert(Geo_Element->NumNodes[j]);
} }
} }
fprintf(PostStream, "$Nodes\n%d\n", (int)nodes.size()); fprintf(PostStream, "$Nodes\n%d\n", (int)nodes.size());
for (int i = 0 ; i < Geo_GetNbrGeoNodes() ; i++) { for(int i = 0; i < Geo_GetNbrGeoNodes(); i++) {
Geo_Node *Geo_Node = Geo_GetGeoNode(i); Geo_Node *Geo_Node = Geo_GetGeoNode(i);
if(nodes.find(Geo_Node->Num) != nodes.end()){ if(nodes.find(Geo_Node->Num) != nodes.end()) {
if(Flag_BIN){ if(Flag_BIN) {
fwrite(&Geo_Node->Num, sizeof(int), 1, PostStream); fwrite(&Geo_Node->Num, sizeof(int), 1, PostStream);
double data[3] = {Geo_Node->x, Geo_Node->y, Geo_Node->z}; double data[3] = {Geo_Node->x, Geo_Node->y, Geo_Node->z};
fwrite(data, sizeof(double), 3, PostStream); fwrite(data, sizeof(double), 3, PostStream);
} }
else{ else {
fprintf(PostStream, "%d %.16g %.16g %.16g\n", fprintf(PostStream, "%d %.16g %.16g %.16g\n", Geo_Node->Num,
Geo_Node->Num, Geo_Node->x, Geo_Node->y, Geo_Node->z) ; Geo_Node->x, Geo_Node->y, Geo_Node->z);
} }
} }
} }
fprintf(PostStream, "$EndNodes\n$Elements\n%d\n", (int)elements.size()); fprintf(PostStream, "$EndNodes\n$Elements\n%d\n", (int)elements.size());
for (unsigned int i = 0 ; i < elements.size() ; i++) { for(unsigned int i = 0; i < elements.size(); i++) {
Geo_Element *Geo_Element = elements[i]; Geo_Element *Geo_Element = elements[i];
int Type = GetDP2Gmsh(Geo_Element->Type) ; int Type = GetDP2Gmsh(Geo_Element->Type);
if(Flag_BIN){ if(Flag_BIN) {
int blob[6] = {Type, 1, 2, Geo_Element->Num, Geo_Element->Region, int blob[6] = {Type,
1,
2,
Geo_Element->Num,
Geo_Element->Region,
Geo_Element->ElementaryRegion}; Geo_Element->ElementaryRegion};
fwrite(blob, sizeof(int), 6, PostStream); fwrite(blob, sizeof(int), 6, PostStream);
std::vector<int> verts(Geo_Element->NbrNodes); std::vector<int> verts(Geo_Element->NbrNodes);
for (int j = 0 ; j < Geo_Element->NbrNodes ; j++) for(int j = 0; j < Geo_Element->NbrNodes; j++)
verts[j] = Geo_Element->NumNodes[j] ; verts[j] = Geo_Element->NumNodes[j];
fwrite(&verts[0], sizeof(int), Geo_Element->NbrNodes, PostStream); fwrite(&verts[0], sizeof(int), Geo_Element->NbrNodes, PostStream);
} }
else{ else {
fprintf(PostStream, "%d %d 2 %d %d ", Geo_Element->Num, fprintf(PostStream, "%d %d 2 %d %d ", Geo_Element->Num, Type,
Type, Geo_Element->Region, Geo_Element->ElementaryRegion) ; Geo_Element->Region, Geo_Element->ElementaryRegion);
for (int j = 0 ; j < Geo_Element->NbrNodes ; j++) for(int j = 0; j < Geo_Element->NbrNodes; j++)
fprintf(PostStream, "%d ", Geo_Element->NumNodes[j]) ; fprintf(PostStream, "%d ", Geo_Element->NumNodes[j]);
fprintf(PostStream, "\n") ; fprintf(PostStream, "\n");
} }
} }
fprintf(PostStream, "$EndElements\n"); fprintf(PostStream, "$EndElements\n");
} }
} }
else if(PostStream && Flag_BIN){/* bricolage */ else if(PostStream && Flag_BIN) { /* bricolage */
fprintf(PostStream, "$PostFormat /* Gmsh 1.2, %s */\n", fprintf(PostStream, "$PostFormat /* Gmsh 1.2, %s */\n",
Flag_BIN ? "binary" : "ascii") ; Flag_BIN ? "binary" : "ascii");
fprintf(PostStream, "1.2 %d %d\n", Flag_BIN, (int)sizeof(double)) ; fprintf(PostStream, "1.2 %d %d\n", Flag_BIN, (int)sizeof(double));
fprintf(PostStream, "$EndPostFormat\n") ; fprintf(PostStream, "$EndPostFormat\n");
} }
break ; break;
case FORMAT_GNUPLOT : case FORMAT_GNUPLOT:
fprintf(PostStream, "# GetDP %s, %s\n", GETDP_VERSION, fprintf(PostStream, "# GetDP %s, %s\n", GETDP_VERSION,
Flag_BIN ? "binary" : "ascii") ; Flag_BIN ? "binary" : "ascii");
break ; break;
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* F o r m a t _ P o s t H e a d e r */ /* F o r m a t _ P o s t H e a d e r */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Format_PostHeader(struct PostSubOperation *PSO_P, int NbTimeStep, void Format_PostHeader(struct PostSubOperation *PSO_P, int NbTimeStep,
int Order, char *Name1, char *Name2) int Order, char *Name1, char *Name2)
{ {
int Format = PSO_P->Format; int Format = PSO_P->Format;
double Time = Current.Time; double Time = Current.Time;
int TimeStep = Current.TimeStep; int TimeStep = Current.TimeStep;
int Contour = PSO_P->Iso; int Contour = PSO_P->Iso;
int Type = PSO_P->CombinationType; int Type = PSO_P->CombinationType;
char name[256] ; char name[256];
CurrentPartitionNumber = 0; CurrentPartitionNumber = 0;
if(Contour){ if(Contour) {
if(!PostElement_L) if(!PostElement_L)
PostElement_L = List_Create(20, 20, sizeof(struct PostElement*)); PostElement_L = List_Create(20, 20, sizeof(struct PostElement *));
else else
List_Reset(PostElement_L); List_Reset(PostElement_L);
} }
if(Name1 && Name2) { if(Name1 && Name2) {
strcpy(name, Order ? Name1 : Name2) ; strcpy(name, Order ? Name1 : Name2);
strcat(name, Get_StringForDefine(PostSubOperation_CombinationType, Type)) ; strcat(name, Get_StringForDefine(PostSubOperation_CombinationType, Type));
strcat(name, Order ? Name2 : Name1) ; strcat(name, Order ? Name2 : Name1);
} }
else if(Name1) else if(Name1)
strcpy(name, Name1) ; strcpy(name, Name1);
else if(Name2) else if(Name2)
strcpy(name, Name2) ; strcpy(name, Name2);
else else
strcpy(name, "unnamed"); strcpy(name, "unnamed");
strcpy(CurrentName, name); strcpy(CurrentName, name);
switch(Format){ switch(Format) {
case FORMAT_GMSH_PARSED : case FORMAT_GMSH_PARSED:
if(PostStream) fprintf(PostStream, "View \"%s\" {\n", name) ; if(PostStream) fprintf(PostStream, "View \"%s\" {\n", name);
Gmsh_StartNewView = 1 ; Gmsh_StartNewView = 1;
break ; break;
case FORMAT_GMSH : case FORMAT_GMSH:
Gmsh_StartNewView = 1 ; Gmsh_StartNewView = 1;
if((PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0) && if((PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0) &&
!PSO_P->FileOut) break; !PSO_P->FileOut)
if(PostStream && Flag_GMSH_VERSION != 2){ break;
if(Flag_BIN){ /* bricolage */ if(PostStream && Flag_GMSH_VERSION != 2) {
if(Flag_BIN) { /* bricolage */
fprintf(PostStream, "$View /* %s */\n", name); fprintf(PostStream, "$View /* %s */\n", name);
fprintf(PostStream, "%s ", name); fprintf(PostStream, "%s ", name);
} }
else { else {
fprintf(PostStream, "View \"%s\" {\n", name) ; fprintf(PostStream, "View \"%s\" {\n", name);
} }
} }
break ; break;
case FORMAT_NXUNV : case FORMAT_NXUNV:
if(PostStream){ if(PostStream) {
NodeTable_StartNew = 1 ; NodeTable_StartNew = 1;
Unv_PrintHeader(PostStream, name, Time, TimeStep, NXUnv_UnitFac Unv_PrintHeader(PostStream, name, Time, TimeStep, NXUnv_UnitFactor,
tor, NXUnv_DatasetLocation); NXUnv_DatasetLocation);
} }
break ; break;
case FORMAT_GNUPLOT : case FORMAT_GNUPLOT:
if(PostStream){ if(PostStream) {
fprintf(PostStream, "# PostData '%s'\n", name); fprintf(PostStream, "# PostData '%s'\n", name);
fprintf(PostStream, "# Type Num X Y Z N1 N2 N3 Values <Values>...\n"); fprintf(PostStream, "# Type Num X Y Z N1 N2 N3 Values <Values>...\n");
} }
break ; break;
case FORMAT_NODE_TABLE : case FORMAT_NODE_TABLE: NodeTable_StartNew = 1; break;
NodeTable_StartNew = 1 ; case FORMAT_ELEMENT_TABLE: ElementTable_StartNew = 1; break;
break ; case FORMAT_SPACE_TABLE:
case FORMAT_ELEMENT_TABLE : case FORMAT_TIME_TABLE:
ElementTable_StartNew = 1 ; case FORMAT_SIMPLE_SPACE_TABLE:
break ; case FORMAT_VALUE_ONLY: TableList_StartNew = 1; break;
case FORMAT_SPACE_TABLE : case FORMAT_ADAPT:
case FORMAT_TIME_TABLE : if(PostStream) fprintf(PostStream, "$Adapt /* %s */\n", name);
case FORMAT_SIMPLE_SPACE_TABLE : break;
case FORMAT_VALUE_ONLY :
TableList_StartNew = 1 ;
break ;
case FORMAT_ADAPT :
if(PostStream) fprintf(PostStream, "$Adapt /* %s */\n", name) ;
break ;
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* F o r m a t _ P o s t F o o t e r */ /* F o r m a t _ P o s t F o o t e r */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Format_PostFooter(struct PostSubOperation *PSO_P, int Store, void Format_PostFooter(struct PostSubOperation *PSO_P, int Store,
bool SendToServer) bool SendToServer)
{ {
List_T *Iso_L[NBR_MAX_ISO], *Solutions_L; List_T *Iso_L[NBR_MAX_ISO], *Solutions_L;
double IsoMin = 1.e200, IsoMax = -1.e200, IsoVal = 0.0, freq, valr, vali ; double IsoMin = 1.e200, IsoMax = -1.e200, IsoVal = 0.0, freq, valr, vali;
int NbrIso = 0 ; int NbrIso = 0;
int iPost, iNode, iIso, iTime, One=1, i, j, NbTimeStep ; int iPost, iNode, iIso, iTime, One = 1, i, j, NbTimeStep;
char tmp[1024]; char tmp[1024];
bool PostOpSolutionGenerated; bool PostOpSolutionGenerated;
struct PostElement *PE ; struct PostElement *PE;
struct Solution *Solution_P=NULL, Solution_S; struct Solution *Solution_P = NULL, Solution_S;
if( !(NbTimeStep = List_Nbr(PSO_P->TimeStep_L)) ) if(!(NbTimeStep = List_Nbr(PSO_P->TimeStep_L)))
NbTimeStep = List_Nbr(Current.DofData->Solutions); NbTimeStep = List_Nbr(Current.DofData->Solutions);
if ( (PSO_P->Format == FORMAT_GMSH || PSO_P->Format == FORMAT_GMSH_PARSED) && if((PSO_P->Format == FORMAT_GMSH || PSO_P->Format == FORMAT_GMSH_PARSED) &&
Flag_GMSH_VERSION != 2 ){ Flag_GMSH_VERSION != 2) {
switch(PSO_P->Legend) {
switch(PSO_P->Legend){
case LEGEND_TIME: case LEGEND_TIME:
Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0], Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0],
PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]); PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]);
for (i = 0 ; i < NbTimeStep ; i++) { for(i = 0; i < NbTimeStep; i++) {
Pos_InitAllSolutions(PSO_P->TimeStep_L, i) ; Pos_InitAllSolutions(PSO_P->TimeStep_L, i);
valr = Current.DofData->CurrentSolution->Time ; valr = Current.DofData->CurrentSolution->Time;
for (j = 0 ; j < Current.NbrHar ; j++){ for(j = 0; j < Current.NbrHar; j++) {
sprintf(tmp, "Step %d/%d: Time = %g", i+1, NbTimeStep, valr); sprintf(tmp, "Step %d/%d: Time = %g", i + 1, NbTimeStep, valr);
Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp); Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp);
} }
} }
Gmsh_StringEnd(PSO_P->Format); Gmsh_StringEnd(PSO_P->Format);
break; break;
case LEGEND_FREQUENCY: case LEGEND_FREQUENCY:
if(Current.NbrHar > 1) { if(Current.NbrHar > 1) {
Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0], Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0],
PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]); PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]);
for (i = 0 ; i < NbTimeStep ; i++) { for(i = 0; i < NbTimeStep; i++) {
Pos_InitAllSolutions(PSO_P->TimeStep_L, i) ; Pos_InitAllSolutions(PSO_P->TimeStep_L, i);
for (j = 0 ; j < Current.NbrHar ; j+=2) { for(j = 0; j < Current.NbrHar; j += 2) {
freq = 0.5*Current.DofData->Val_Pulsation[j/2]/M_PI ; freq = 0.5 * Current.DofData->Val_Pulsation[j / 2] / M_PI;
sprintf(tmp, "%g Hz (real part: COSINUS)", freq); sprintf(tmp, "%g Hz (real part: COSINUS)", freq);
Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp); Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp);
sprintf(tmp, "%g Hz (imag part: -SINUS)", freq); sprintf(tmp, "%g Hz (imag part: -SINUS)", freq);
Gmsh_StringAdd(PSO_P->Format, 0, tmp); Gmsh_StringAdd(PSO_P->Format, 0, tmp);
} }
} }
Gmsh_StringEnd(PSO_P->Format); Gmsh_StringEnd(PSO_P->Format);
} }
break; break;
case LEGEND_EIGENVALUES: case LEGEND_EIGENVALUES:
Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0], Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0],
PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]); PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]);
for (i = 0 ; i < NbTimeStep ; i++) { for(i = 0; i < NbTimeStep; i++) {
Pos_InitAllSolutions(PSO_P->TimeStep_L, i) ; Pos_InitAllSolutions(PSO_P->TimeStep_L, i);
valr = Current.DofData->CurrentSolution->Time ; valr = Current.DofData->CurrentSolution->Time;
vali = Current.DofData->CurrentSolution->TimeImag ; vali = Current.DofData->CurrentSolution->TimeImag;
if(Current.NbrHar == 1){ if(Current.NbrHar == 1) {
sprintf(tmp, "Eigenvalue %d/%d: %g", sprintf(tmp, "Eigenvalue %d/%d: %g", i + 1, NbTimeStep, valr);
i+1, NbTimeStep, valr);
Gmsh_StringAdd(PSO_P->Format, !i, tmp); Gmsh_StringAdd(PSO_P->Format, !i, tmp);
} }
else{ else {
for (j = 0 ; j < Current.NbrHar ; j++) { for(j = 0; j < Current.NbrHar; j++) {
if(!(j % 2)){ if(!(j % 2)) {
sprintf(tmp, "Eigenvalue %d/%d: %g %s i * %g (Real Part)", sprintf(tmp, "Eigenvalue %d/%d: %g %s i * %g (Real Part)", i + 1,
i+1, NbTimeStep, valr, (vali > 0) ? "+" : "-", NbTimeStep, valr, (vali > 0) ? "+" : "-",
(vali > 0) ? vali : -vali); (vali > 0) ? vali : -vali);
Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp); Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp);
} }
else{ else {
sprintf(tmp, "Eigenvalue %d/%d: %g %s i * %g (Imaginary Part)", sprintf(tmp, "Eigenvalue %d/%d: %g %s i * %g (Imaginary Part)",
i+1, NbTimeStep, valr, (vali > 0) ? "+" : "-", i + 1, NbTimeStep, valr, (vali > 0) ? "+" : "-",
(vali > 0) ? vali : -vali); (vali > 0) ? vali : -vali);
Gmsh_StringAdd(PSO_P->Format, 0, tmp); Gmsh_StringAdd(PSO_P->Format, 0, tmp);
} }
} }
} }
} }
Gmsh_StringEnd(PSO_P->Format); Gmsh_StringEnd(PSO_P->Format);
break; break;
} }
} }
if(PSO_P->Iso){ if(PSO_P->Iso) {
for(iPost = 0; iPost < List_Nbr(PostElement_L); iPost++) {
for(iPost = 0 ; iPost < List_Nbr(PostElement_L) ; iPost++){ PE = *(struct PostElement **)List_Pointer(PostElement_L, iPost);
PE = *(struct PostElement**)List_Pointer(PostElement_L, iPost); for(iNode = 0; iNode < PE->NbrNodes; iNode++) {
for (iNode = 0 ; iNode < PE->NbrNodes ; iNode++ ){ IsoMin = std::min(IsoMin, PE->Value[iNode].Val[0]);
IsoMin = std::min(IsoMin, PE->Value[iNode].Val[0]) ; IsoMax = std::max(IsoMax, PE->Value[iNode].Val[0]);
IsoMax = std::max(IsoMax, PE->Value[iNode].Val[0]) ;
} }
} }
if((NbrIso = PSO_P->Iso) < 0) if((NbrIso = PSO_P->Iso) < 0) NbrIso = List_Nbr(PSO_P->Iso_L);
NbrIso = List_Nbr(PSO_P->Iso_L) ; if(NbrIso > NBR_MAX_ISO) {
if(NbrIso > NBR_MAX_ISO){
Message::Error("Too many Iso values"); Message::Error("Too many Iso values");
NbrIso = NBR_MAX_ISO; NbrIso = NBR_MAX_ISO;
} }
if(PostStream && PSO_P->Format == FORMAT_GNUPLOT) if(PostStream && PSO_P->Format == FORMAT_GNUPLOT)
fprintf(PostStream, "# NbIso = %d, Min = %g, Max = %g\n", fprintf(PostStream, "# NbIso = %d, Min = %g, Max = %g\n", NbrIso, IsoMin,
NbrIso, IsoMin, IsoMax) ; IsoMax);
for(iIso = 0 ; iIso < NbrIso ; iIso++) for(iIso = 0; iIso < NbrIso; iIso++)
Iso_L[iIso] = List_Create(10, 10, sizeof(struct PostElement*)) ; Iso_L[iIso] = List_Create(10, 10, sizeof(struct PostElement *));
for(iPost = 0 ; iPost < List_Nbr(PostElement_L) ; iPost++){ for(iPost = 0; iPost < List_Nbr(PostElement_L); iPost++) {
PE = *(struct PostElement**)List_Pointer(PostElement_L, iPost); PE = *(struct PostElement **)List_Pointer(PostElement_L, iPost);
for(iIso = 0 ; iIso < NbrIso ; iIso++){ for(iIso = 0; iIso < NbrIso; iIso++) {
if(PSO_P->Iso > 0){ if(PSO_P->Iso > 0) {
Cal_Iso(PE, Iso_L[iIso], IsoMin+iIso*(IsoMax-IsoMin)/(double)(NbrIso-1) Cal_Iso(PE, Iso_L[iIso],
, IsoMin + iIso * (IsoMax - IsoMin) / (double)(NbrIso - 1),
IsoMin, IsoMax, PSO_P->DecomposeInSimplex) ; IsoMin, IsoMax, PSO_P->DecomposeInSimplex);
} }
else{ else {
List_Read(PSO_P->Iso_L, iIso, &IsoVal) ; List_Read(PSO_P->Iso_L, iIso, &IsoVal);
Cal_Iso(PE, Iso_L[iIso], IsoVal, IsoMin, IsoMax, PSO_P->DecomposeInSimp Cal_Iso(PE, Iso_L[iIso], IsoVal, IsoMin, IsoMax,
lex) ; PSO_P->DecomposeInSimplex);
} }
} }
if(!Store) Destroy_PostElement(PE); if(!Store) Destroy_PostElement(PE);
} }
for(iIso = 0 ; iIso < NbrIso ; iIso++){ for(iIso = 0; iIso < NbrIso; iIso++) {
for(iPost = 0 ; iPost < List_Nbr(Iso_L[iIso]) ; iPost++){ for(iPost = 0; iPost < List_Nbr(Iso_L[iIso]); iPost++) {
PE = *(struct PostElement**)List_Pointer(Iso_L[iIso], iPost) ; PE = *(struct PostElement **)List_Pointer(Iso_L[iIso], iPost);
Format_PostElement(PSO_P, 0, 0, Format_PostElement(PSO_P, 0, 0, Current.Time, 0, 1, Current.NbrHar,
Current.Time, 0, 1, PSO_P->HarmonicToTime, NULL, PE);
Current.NbrHar, PSO_P->HarmonicToTime, Destroy_PostElement(PE);
NULL, PE); }
Destroy_PostElement(PE) ; List_Delete(Iso_L[iIso]);
} if(PostStream && PSO_P->Format == FORMAT_GNUPLOT)
List_Delete(Iso_L[iIso]) ; fprintf(PostStream, "\n");
if(PostStream && PSO_P->Format == FORMAT_GNUPLOT) fprintf(PostStream, "\n"
) ;
} }
} }
switch(PSO_P->Format){ switch(PSO_P->Format) {
case FORMAT_GMSH_PARSED : case FORMAT_GMSH_PARSED:
if(PostStream && List_Nbr(TimeValue_L) > 1){ if(PostStream && List_Nbr(TimeValue_L) > 1) {
fprintf(PostStream, "TIME{"); fprintf(PostStream, "TIME{");
for(iTime = 0; iTime < List_Nbr(TimeValue_L); iTime++){ for(iTime = 0; iTime < List_Nbr(TimeValue_L); iTime++) {
if(iTime) fprintf(PostStream, ","); if(iTime) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g", *(double*)List_Pointer(TimeValue_L, iTime)); fprintf(PostStream, "%.16g",
*(double *)List_Pointer(TimeValue_L, iTime));
} }
fprintf(PostStream, "};\n"); fprintf(PostStream, "};\n");
} }
fprintf(PostStream, "};\n") ; fprintf(PostStream, "};\n");
break ; break;
case FORMAT_GMSH : case FORMAT_GMSH:
if(Gmsh_StartNewView) Gmsh_ResetStaticLists(); // nothing to print! if(Gmsh_StartNewView) Gmsh_ResetStaticLists(); // nothing to print!
if(PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0){ if(PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0) {
#if defined(HAVE_GMSH) #if defined(HAVE_GMSH)
int field = (PSO_P->StoreInField >= 0) ? PSO_P->StoreInField : int field = (PSO_P->StoreInField >= 0) ? PSO_P->StoreInField :
PSO_P->StoreInMeshBasedField; PSO_P->StoreInMeshBasedField;
Message::Info("Storing data in field %d (%s)", field, Message::Info("Storing data in field %d (%s)", field,
PSO_P->StoreInField >= 0 ? "list-based" : "mesh-based"); PSO_P->StoreInField >= 0 ? "list-based" : "mesh-based");
int NS[24] = {NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT, NbTT, int NS[24] = {NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT,
NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH, NbTT, NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH,
NbSI, NbVI, NbTI, NbSY, NbVY, NbTY}; NbVH, NbTH, NbSI, NbVI, NbTI, NbSY, NbVY, NbTY};
std::vector<double> *LS[24] = {&SP, &VP, &TP, &SL, &VL, &TL, &ST, &VT, & std::vector<double> *LS[24] = {&SP, &VP, &TP, &SL, &VL, &TL, &ST, &VT,
TT, &TT, &SQ, &VQ, &TQ, &SS, &VS, &TS1, &SH,
&SQ, &VQ, &TQ, &SS, &VS, &TS1, &SH, &VH, &VH, &TH, &SI, &VI, &TI, &SY, &VY, &TY};
&TH,
&SI, &VI, &TI, &SY, &VY, &TY};
PViewData *data; PViewData *data;
if(PSO_P->StoreInField >= 0) if(PSO_P->StoreInField >= 0)
data = new PViewDataList(); data = new PViewDataList();
else else
data = new PViewDataGModel(PViewDataGModel::ElementNodeData); data = new PViewDataGModel(PViewDataGModel::ElementNodeData);
data->importLists(NS, LS); data->importLists(NS, LS);
new PView(data, field); new PView(data, field);
#else #else
Message::Error("GetDP must be compiled with Gmsh support to store data as Message::Error(
field"); "GetDP must be compiled with Gmsh support to store data as field");
#endif #endif
if(!PSO_P->FileOut) break; if(!PSO_P->FileOut) break;
} }
if(Flag_GMSH_VERSION == 2){ if(Flag_GMSH_VERSION == 2) {
int NS[8] = {NbSP, NbSL, NbST, NbSQ, NbSS, NbSH, NbSI, NbSY}; int NS[8] = {NbSP, NbSL, NbST, NbSQ, NbSS, NbSH, NbSI, NbSY};
std::vector<double> *LS[8] = {&SP, &SL, &ST, &SQ, &SS, &SH, &SI, &SY}; std::vector<double> *LS[8] = {&SP, &SL, &ST, &SQ, &SS, &SH, &SI, &SY};
Gmsh_PrintElementNodeData(PSO_P, NbTimeStep, 1, NS, LS); Gmsh_PrintElementNodeData(PSO_P, NbTimeStep, 1, NS, LS);
int NV[8] = {NbVP, NbVL, NbVT, NbVQ, NbVS, NbVH, NbVI, NbVY}; int NV[8] = {NbVP, NbVL, NbVT, NbVQ, NbVS, NbVH, NbVI, NbVY};
std::vector<double> *LV[8] = {&VP, &VL, &VT, &VQ, &VS, &VH, &VI, &VY}; std::vector<double> *LV[8] = {&VP, &VL, &VT, &VQ, &VS, &VH, &VI, &VY};
Gmsh_PrintElementNodeData(PSO_P, NbTimeStep, 3, NV, LV); Gmsh_PrintElementNodeData(PSO_P, NbTimeStep, 3, NV, LV);
int NT[8] = {NbTP, NbTL, NbTT, NbTQ, NbTS, NbTH, NbTI, NbTY}; int NT[8] = {NbTP, NbTL, NbTT, NbTQ, NbTS, NbTH, NbTI, NbTY};
std::vector<double> *LT[8] = {&TP, &TL, &TT, &TQ, &TS1, &TH, &TI, &TY}; std::vector<double> *LT[8] = {&TP, &TL, &TT, &TQ, &TS1, &TH, &TI, &TY};
Gmsh_PrintElementNodeData(PSO_P, NbTimeStep, 9, NT, LT); Gmsh_PrintElementNodeData(PSO_P, NbTimeStep, 9, NT, LT);
} }
else if(PostStream && Flag_BIN){ /* bricolage */ else if(PostStream && Flag_BIN) { /* bricolage */
fprintf(PostStream, "%d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d " fprintf(PostStream,
"%d %d %d %d %d %d %d %d %d %d %d 0 0\n", "%d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d "
List_Nbr(TimeValue_L), "%d %d %d %d %d %d %d %d %d %d %d 0 0\n",
NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT, NbTT, List_Nbr(TimeValue_L), NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST,
NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH, NbVT, NbTT, NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH,
NbSI, NbVI, NbTI, NbSY, NbVY, NbTY, NbT2, (int)T2C.size()); NbSI, NbVI, NbTI, NbSY, NbVY, NbTY, NbT2, (int)T2C.size());
fwrite(&One, sizeof(int), 1, PostStream); fwrite(&One, sizeof(int), 1, PostStream);
List_WriteToFile(TimeValue_L, PostStream, LIST_FORMAT_BINARY); List_WriteToFile(TimeValue_L, PostStream, LIST_FORMAT_BINARY);
bool f = true; bool f = true;
dVecWrite(SP, PostStream, f); dVecWrite(VP, PostStream, f); dVecWrite(SP, PostStream, f);
dVecWrite(VP, PostStream, f);
dVecWrite(TP, PostStream, f); dVecWrite(TP, PostStream, f);
dVecWrite(SL, PostStream, f); dVecWrite(VL, PostStream, f); dVecWrite(SL, PostStream, f);
dVecWrite(VL, PostStream, f);
dVecWrite(TL, PostStream, f); dVecWrite(TL, PostStream, f);
dVecWrite(ST, PostStream, f); dVecWrite(VT, PostStream, f); dVecWrite(ST, PostStream, f);
dVecWrite(VT, PostStream, f);
dVecWrite(TT, PostStream, f); dVecWrite(TT, PostStream, f);
dVecWrite(SQ, PostStream, f); dVecWrite(VQ, PostStream, f); dVecWrite(SQ, PostStream, f);
dVecWrite(VQ, PostStream, f);
dVecWrite(TQ, PostStream, f); dVecWrite(TQ, PostStream, f);
dVecWrite(SS, PostStream, f); dVecWrite(VS, PostStream, f); dVecWrite(SS, PostStream, f);
dVecWrite(VS, PostStream, f);
dVecWrite(TS1, PostStream, f); dVecWrite(TS1, PostStream, f);
dVecWrite(SH, PostStream, f); dVecWrite(VH, PostStream, f); dVecWrite(SH, PostStream, f);
dVecWrite(VH, PostStream, f);
dVecWrite(TH, PostStream, f); dVecWrite(TH, PostStream, f);
dVecWrite(SI, PostStream, f); dVecWrite(VI, PostStream, f); dVecWrite(SI, PostStream, f);
dVecWrite(VI, PostStream, f);
dVecWrite(TI, PostStream, f); dVecWrite(TI, PostStream, f);
dVecWrite(SY, PostStream, f); dVecWrite(VY, PostStream, f); dVecWrite(SY, PostStream, f);
dVecWrite(VY, PostStream, f);
dVecWrite(TY, PostStream, f); dVecWrite(TY, PostStream, f);
dVecWrite(T2D, PostStream, f); cVecWrite(T2C, PostStream, f); dVecWrite(T2D, PostStream, f);
cVecWrite(T2C, PostStream, f);
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
fprintf(PostStream, "$EndView\n"); fprintf(PostStream, "$EndView\n");
} }
else if(PostStream){ else if(PostStream) {
if(List_Nbr(TimeValue_L) > 1){ if(List_Nbr(TimeValue_L) > 1) {
fprintf(PostStream, "TIME{"); fprintf(PostStream, "TIME{");
for(iTime = 0; iTime < List_Nbr(TimeValue_L); iTime++){ for(iTime = 0; iTime < List_Nbr(TimeValue_L); iTime++) {
if(iTime) fprintf(PostStream, ","); if(iTime) fprintf(PostStream, ",");
fprintf(PostStream, "%.16g", *(double*)List_Pointer(TimeValue_L, iTime) fprintf(PostStream, "%.16g",
); *(double *)List_Pointer(TimeValue_L, iTime));
} }
fprintf(PostStream, "};\n"); fprintf(PostStream, "};\n");
} }
fprintf(PostStream, "};\n") ; fprintf(PostStream, "};\n");
} }
break ; break;
case FORMAT_ADAPT : case FORMAT_ADAPT:
if(PostStream) fprintf(PostStream, "$EndAdapt\n"); if(PostStream) fprintf(PostStream, "$EndAdapt\n");
break ; break;
case FORMAT_NXUNV : case FORMAT_NXUNV:
if(PostStream){ if(PostStream) {
if(NXUnv_DatasetLocation == 1) //Data at nodes if(NXUnv_DatasetLocation == 1) // Data at nodes
Unv_PrintNodeTable(PostStream, NodeTable, NXUnv_UnitFac Unv_PrintNodeTable(PostStream, NodeTable, NXUnv_UnitFactor);
tor); Unv_PrintFooter(PostStream);
Unv_PrintFooter(PostStream); }
} break;
break ; case FORMAT_NODE_TABLE:
case FORMAT_NODE_TABLE : if(PostStream && NodeTable.size()) {
if(PostStream && NodeTable.size()){
fprintf(PostStream, "%d\n", (int)NodeTable.size()); fprintf(PostStream, "%d\n", (int)NodeTable.size());
for(std::map<int, std::vector<double> >::iterator it = NodeTable.begin(); for(std::map<int, std::vector<double> >::iterator it = NodeTable.begin();
it != NodeTable.end(); it++){ it != NodeTable.end(); it++) {
fprintf(PostStream, "%d", it->first); fprintf(PostStream, "%d", it->first);
for(unsigned int i = 0; i < it->second.size(); i++) for(unsigned int i = 0; i < it->second.size(); i++)
fprintf(PostStream, " %.16g", it->second[i]); fprintf(PostStream, " %.16g", it->second[i]);
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
} }
} }
{ {
std::vector<double> exp; std::vector<double> exp;
exp.push_back(NodeTable.size()); exp.push_back(NodeTable.size());
for(std::map<int, std::vector<double> >::iterator it = NodeTable.begin(); for(std::map<int, std::vector<double> >::iterator it = NodeTable.begin();
it != NodeTable.end(); it++){ it != NodeTable.end(); it++) {
exp.push_back(it->first); exp.push_back(it->first);
for(unsigned int i = 0; i < it->second.size(); i++) for(unsigned int i = 0; i < it->second.size(); i++)
exp.push_back(it->second[i]); exp.push_back(it->second[i]);
} }
GetDPNumbers[CurrentName] = exp; GetDPNumbers[CurrentName] = exp;
GetDPNumbersMap[CurrentName] = NodeTable; GetDPNumbersMap[CurrentName] = NodeTable;
if(SendToServer && PSO_P->SendToServer && strcmp(PSO_P->SendToServer, "No" if(SendToServer && PSO_P->SendToServer &&
)) strcmp(PSO_P->SendToServer, "No"))
Message::AddOnelabNumberChoice(PSO_P->SendToServer, exp, PSO_P->Color, Message::AddOnelabNumberChoice(PSO_P->SendToServer, exp, PSO_P->Color,
PSO_P->Units, PSO_P->Label, PSO_P->Visibl PSO_P->Units, PSO_P->Label,
e, PSO_P->Visible, PSO_P->Closed);
PSO_P->Closed);
} }
break; break;
case FORMAT_ELEMENT_TABLE : case FORMAT_ELEMENT_TABLE:
if(PostStream){ if(PostStream) {
fprintf(PostStream, "%d\n", (int)ElementTable.size()); fprintf(PostStream, "%d\n", (int)ElementTable.size());
for(std::map<int, std::vector<double> >::iterator it = ElementTable.begin( for(std::map<int, std::vector<double> >::iterator it =
); ElementTable.begin();
it != ElementTable.end(); it++){ it != ElementTable.end(); it++) {
fprintf(PostStream, "%d", it->first); fprintf(PostStream, "%d", it->first);
for(unsigned int i = 0; i < it->second.size(); i++) for(unsigned int i = 0; i < it->second.size(); i++)
fprintf(PostStream, " %.16g", it->second[i]); fprintf(PostStream, " %.16g", it->second[i]);
fprintf(PostStream, "\n"); fprintf(PostStream, "\n");
} }
} }
{ {
std::vector<double> exp; std::vector<double> exp;
exp.push_back(ElementTable.size()); exp.push_back(ElementTable.size());
for(std::map<int, std::vector<double> >::iterator it = ElementTable.begin( for(std::map<int, std::vector<double> >::iterator it =
); ElementTable.begin();
it != ElementTable.end(); it++){ it != ElementTable.end(); it++) {
exp.push_back(it->first); exp.push_back(it->first);
for(unsigned int i = 0; i < it->second.size(); i++) for(unsigned int i = 0; i < it->second.size(); i++)
exp.push_back(it->second[i]); exp.push_back(it->second[i]);
} }
GetDPNumbers[CurrentName] = exp; GetDPNumbers[CurrentName] = exp;
GetDPNumbersMap[CurrentName] = ElementTable; GetDPNumbersMap[CurrentName] = ElementTable;
if(SendToServer && PSO_P->SendToServer && strcmp(PSO_P->SendToServer, "No" if(SendToServer && PSO_P->SendToServer &&
)) strcmp(PSO_P->SendToServer, "No"))
Message::AddOnelabNumberChoice(PSO_P->SendToServer, exp, PSO_P->Color, Message::AddOnelabNumberChoice(PSO_P->SendToServer, exp, PSO_P->Color,
PSO_P->Units, PSO_P->Label, PSO_P->Visibl PSO_P->Units, PSO_P->Label,
e, PSO_P->Visible, PSO_P->Closed);
PSO_P->Closed);
}
break;
case FORMAT_SPACE_TABLE :
case FORMAT_TIME_TABLE :
case FORMAT_SIMPLE_SPACE_TABLE :
case FORMAT_VALUE_ONLY :
{
if(TableList.size()) {
std::vector<double> v(TableList.begin(), TableList.end());
GetDPNumbers[CurrentName] = v;
if(SendToServer && PSO_P->SendToServer && strcmp(PSO_P->SendToServer, "N
o"))
Message::AddOnelabNumberChoice(PSO_P->SendToServer, v, PSO_P->Color,
PSO_P->Units, PSO_P->Label, PSO_P->Visi
ble,
PSO_P->Closed);
}
} }
break; break;
case FORMAT_LOOP_ERROR : case FORMAT_SPACE_TABLE:
Solutions_L = ((struct PostOpSolutions*) case FORMAT_TIME_TABLE:
List_Pointer(Current.PostOpData_L, Current.PostOpDataIndex))->Solutions_L; case FORMAT_SIMPLE_SPACE_TABLE:
case FORMAT_VALUE_ONLY: {
if(TableList.size()) {
std::vector<double> v(TableList.begin(), TableList.end());
GetDPNumbers[CurrentName] = v;
if(SendToServer && PSO_P->SendToServer &&
strcmp(PSO_P->SendToServer, "No"))
Message::AddOnelabNumberChoice(PSO_P->SendToServer, v, PSO_P->Color,
PSO_P->Units, PSO_P->Label,
PSO_P->Visible, PSO_P->Closed);
}
} break;
case FORMAT_LOOP_ERROR:
Solutions_L = ((struct PostOpSolutions *)List_Pointer(
Current.PostOpData_L, Current.PostOpDataIndex))
->Solutions_L;
PostOpSolutionGenerated = false; PostOpSolutionGenerated = false;
if(List_Nbr(Solutions_L)>0) { if(List_Nbr(Solutions_L) > 0) {
Solution_P = (struct Solution*)List_Pointer(Solutions_L, List_Nbr(Solution Solution_P =
s_L)-1); (struct Solution *)List_Pointer(Solutions_L, List_Nbr(Solutions_L) - 1);
PostOpSolutionGenerated = (Solution_P->TimeStep == (int)Current.TimeStep); PostOpSolutionGenerated = (Solution_P->TimeStep == (int)Current.TimeStep);
} }
if (!PostOpSolutionGenerated) { if(!PostOpSolutionGenerated) {
Solution_S.Time = Current.Time; Solution_S.Time = Current.Time;
Solution_S.TimeImag = Current.TimeImag; Solution_S.TimeImag = Current.TimeImag;
Solution_S.TimeStep = Current.TimeStep; Solution_S.TimeStep = Current.TimeStep;
Solution_S.SolutionExist = 1; Solution_S.SolutionExist = 1;
Solution_S.TimeFunctionValues = NULL; Solution_S.TimeFunctionValues = NULL;
LinAlg_CreateVector(&Solution_S.x, &Current.DofData->Solver, LinAlg_CreateVector(&Solution_S.x, &Current.DofData->Solver,
List_Nbr(PostOpResults_L)); List_Nbr(PostOpResults_L));
for(int i=0; i<List_Nbr(PostOpResults_L); i++){ for(int i = 0; i < List_Nbr(PostOpResults_L); i++) {
List_Read(PostOpResults_L, i, &valr); List_Read(PostOpResults_L, i, &valr);
LinAlg_SetDoubleInVector(valr, &Solution_S.x, i); LinAlg_SetDoubleInVector(valr, &Solution_S.x, i);
} }
List_Add(Solutions_L, &Solution_S); List_Add(Solutions_L, &Solution_S);
} }
else else
for(int i=0; i<List_Nbr(PostOpResults_L); i++){ for(int i = 0; i < List_Nbr(PostOpResults_L); i++) {
List_Read(PostOpResults_L, i, &valr); List_Read(PostOpResults_L, i, &valr);
LinAlg_SetDoubleInVector(valr, &Solution_P->x, i); LinAlg_SetDoubleInVector(valr, &Solution_P->x, i);
} }
List_Delete(PostOpResults_L); List_Delete(PostOpResults_L);
PostOpResults_L = NULL; PostOpResults_L = NULL;
break; break;
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* F o r m a t _ P o s t E l e m e n t */ /* F o r m a t _ P o s t E l e m e n t */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Format_PostElement(struct PostSubOperation *PSO_P, int Contour, int Store, void Format_PostElement(struct PostSubOperation *PSO_P, int Contour, int Store,
double Time, int TimeStep, int NbTimeStep, double Time, int TimeStep, int NbTimeStep,
int NbrHarmonics, int HarmonicToTime, double *Dummy, int NbrHarmonics, int HarmonicToTime, double *Dummy,
struct PostElement * PE) struct PostElement *PE)
{ {
int i, j, k, l, Num_Element ; int i, j, k, l, Num_Element;
struct PostElement * PE2 ; struct PostElement *PE2;
struct Value Value ; struct Value Value;
static int Warning_FirstHarmonic = 0 ; static int Warning_FirstHarmonic = 0;
/* TODO /* TODO
static int Size = 0 ; static int Size = 0 ;
int flag_storeAllTimeResults, indexInTmpValues; int flag_storeAllTimeResults, indexInTmpValues;
static struct Value TmpValue, *TmpValues ; static struct Value TmpValue, *TmpValues ;
static double *Times ; static double *Times ;
struct Value *FourierValues; struct Value *FourierValues;
skipping to change at line 1641 skipping to change at line 1978
case VECTOR : Size = 3 ; break ; case VECTOR : Size = 3 ; break ;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG : Size = 3 ; break ;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM : Size = 6 ; break ;
case TENSOR : Size = 9 ; break ; case TENSOR : Size = 9 ; break ;
default : Size = 9 ; break ; default : Size = 9 ; break ;
} }
} }
*/ */
if(PE->Index != NO_ELEMENT) if(PE->Index != NO_ELEMENT)
Num_Element = Geo_GetGeoElement(PE->Index)->Num ; Num_Element = Geo_GetGeoElement(PE->Index)->Num;
else else
Num_Element = 0 ; Num_Element = 0;
if(Contour){ if(Contour) {
if(PE->Value[0].Type != SCALAR){ if(PE->Value[0].Type != SCALAR) {
Message::Error("Non scalar Element %d in contour creation", Num_Element); Message::Error("Non scalar Element %d in contour creation", Num_Element);
return; return;
} }
if(NbTimeStep != 1){ if(NbTimeStep != 1) {
Message::Error("Contour creation not allowed for multiple time steps"); Message::Error("Contour creation not allowed for multiple time steps");
return; return;
} }
if(Current.NbrHar != 1 && !Warning_FirstHarmonic){ if(Current.NbrHar != 1 && !Warning_FirstHarmonic) {
Message::Warning("Contour creation done only for first harmonic (use Re[] Message::Warning(
or Im[])"); "Contour creation done only for first harmonic (use Re[] or Im[])");
Warning_FirstHarmonic = 1 ; Warning_FirstHarmonic = 1;
} }
if(Store) if(Store)
List_Add(PostElement_L, &PE) ; List_Add(PostElement_L, &PE);
else{ else {
PE2 = PartialCopy_PostElement(PE) ; PE2 = PartialCopy_PostElement(PE);
List_Add(PostElement_L, &PE2) ; List_Add(PostElement_L, &PE2);
} }
return ; return;
} }
if(PSO_P->ChangeOfCoordinates[0] >= 0){ if(PSO_P->ChangeOfCoordinates[0] >= 0) {
for(i=0 ; i<PE->NbrNodes ; i++){ for(i = 0; i < PE->NbrNodes; i++) {
Current.x = PE->x[i]; Current.x = PE->x[i];
Current.y = PE->y[i]; Current.y = PE->y[i];
Current.z = PE->z[i]; Current.z = PE->z[i];
for(j = 0; j<9 ; j++) Current.Val[j] = PE->Value[i].Val[j]; for(j = 0; j < 9; j++) Current.Val[j] = PE->Value[i].Val[j];
Get_ValueOfExpressionByIndex(PSO_P->ChangeOfCoordinates[0], NULL, 0., 0., Get_ValueOfExpressionByIndex(PSO_P->ChangeOfCoordinates[0], NULL, 0., 0.,
0., &Value) ; 0., &Value);
PE->x[i] = Value.Val[0]; PE->x[i] = Value.Val[0];
Get_ValueOfExpressionByIndex(PSO_P->ChangeOfCoordinates[1], NULL, 0., 0., Get_ValueOfExpressionByIndex(PSO_P->ChangeOfCoordinates[1], NULL, 0., 0.,
0., &Value) ; 0., &Value);
PE->y[i] = Value.Val[0]; PE->y[i] = Value.Val[0];
Get_ValueOfExpressionByIndex(PSO_P->ChangeOfCoordinates[2], NULL, 0., 0., Get_ValueOfExpressionByIndex(PSO_P->ChangeOfCoordinates[2], NULL, 0., 0.,
0., &Value) ; 0., &Value);
PE->z[i] = Value.Val[0]; PE->z[i] = Value.Val[0];
} }
} }
if(PSO_P->ChangeOfValues && List_Nbr(PSO_P->ChangeOfValues) > 0){ if(PSO_P->ChangeOfValues && List_Nbr(PSO_P->ChangeOfValues) > 0) {
for(i=0 ; i < PE->NbrNodes ; i++){ for(i = 0; i < PE->NbrNodes; i++) {
Current.x = PE->x[i]; Current.x = PE->x[i];
Current.y = PE->y[i]; Current.y = PE->y[i];
Current.z = PE->z[i]; Current.z = PE->z[i];
for(k=0 ; k<Current.NbrHar ; k++){ for(k = 0; k < Current.NbrHar; k++) {
for(j = 0; j<9 ; j++) Current.Val[j] = PE->Value[i].Val[MAX_DIM*k+j]; for(j = 0; j < 9; j++)
for(l=0 ; l<List_Nbr(PSO_P->ChangeOfValues) ; l++){ Current.Val[j] = PE->Value[i].Val[MAX_DIM * k + j];
Get_ValueOfExpressionByIndex(*(int*)List_Pointer(PSO_P->ChangeOfValues, for(l = 0; l < List_Nbr(PSO_P->ChangeOfValues); l++) {
l), Get_ValueOfExpressionByIndex(
NULL, 0., 0., 0., &Value) ; *(int *)List_Pointer(PSO_P->ChangeOfValues, l), NULL, 0., 0., 0.,
PE->Value[i].Val[MAX_DIM*k+l] = Value.Val[0]; &Value);
} PE->Value[i].Val[MAX_DIM * k + l] = Value.Val[0];
}
} }
} }
} }
switch(PSO_P->Format){ switch(PSO_P->Format) {
case FORMAT_GMSH_PARSED : case FORMAT_GMSH_PARSED:
GmshParsed_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics, HarmonicTo GmshParsed_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics,
Time, HarmonicToTime, PE->Type, PE->NbrNodes, PE->x,
PE->Type, PE->NbrNodes, PE->x, PE->y, PE->z, PE->y, PE->z, PE->Value);
PE->Value) ; break;
break ; case FORMAT_NXUNV:
case FORMAT_NXUNV : if(PostStream) {
if(PostStream){ if(NXUnv_DatasetLocation == 1) // Data at nodes
if(NXUnv_DatasetLocation == 1) //Data at nodes NodeTable_PrintElement(TimeStep, NbTimeStep, NbrHarmonics, PE);
NodeTable_PrintElement(TimeStep, NbTimeStep, NbrHarmoni else if(NXUnv_DatasetLocation == 2 ||
cs, PE); NXUnv_DatasetLocation ==
else if(NXUnv_DatasetLocation == 2 || NXUnv_DatasetLocation == 3) // Data at elements or nodes on elements
3) //Data at elements or nodes on elements Unv_PrintElement(PostStream, Num_Element, PE->NbrNodes, PE->Value,
Unv_PrintElement(PostStream, Num_Element, PE->NbrNodes, NbrHarmonics, NXUnv_DatasetLocation, NXUnv_UnitFactor);
PE->Value, NbrHarmonics, NXUnv_DatasetLocation, NXUnv_UnitFactor) ; }
} break;
break ; case FORMAT_GMSH:
case FORMAT_GMSH : if(PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0) {
if(PSO_P->StoreInField >= 0 || PSO_P->StoreInMeshBasedField >= 0){ Gmsh_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics,
Gmsh_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics, HarmonicToTime HarmonicToTime, PE->Type, Num_Element, PE->NbrNodes,
, PE->x, PE->y, PE->z, PE->Value, PSO_P, Store);
PE->Type, Num_Element, PE->NbrNodes, PE->x, PE->y, PE->z
,
PE->Value, PSO_P, Store) ;
if(!PSO_P->FileOut || Flag_GMSH_VERSION == 2 || Flag_BIN) break; if(!PSO_P->FileOut || Flag_GMSH_VERSION == 2 || Flag_BIN) break;
} }
if(Flag_GMSH_VERSION == 2 || Flag_BIN){ /* bricolage */ if(Flag_GMSH_VERSION == 2 || Flag_BIN) { /* bricolage */
Gmsh_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics, HarmonicToTime Gmsh_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics,
, HarmonicToTime, PE->Type, Num_Element, PE->NbrNodes,
PE->Type, Num_Element, PE->NbrNodes, PE->x, PE->y, PE->z PE->x, PE->y, PE->z, PE->Value, PSO_P, Store);
, }
PE->Value, PSO_P, Store) ; else {
} GmshParsed_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics,
else{ HarmonicToTime, PE->Type, PE->NbrNodes, PE->x,
GmshParsed_PrintElement(Time, TimeStep, NbTimeStep, NbrHarmonics, Harmonic PE->y, PE->z, PE->Value);
ToTime, }
PE->Type, PE->NbrNodes, PE->x, PE->y, PE->z, break;
PE->Value) ; case FORMAT_GNUPLOT:
} Gnuplot_PrintElement(PSO_P->Format, Time, TimeStep, NbTimeStep,
break ; NbrHarmonics, HarmonicToTime, PE->Type, Num_Element,
case FORMAT_GNUPLOT : PE->NbrNodes, PE->x, PE->y, PE->z, Dummy, PE->Value);
Gnuplot_PrintElement(PSO_P->Format, Time, TimeStep, NbTimeStep, NbrHarmonics break;
, case FORMAT_SPACE_TABLE:
HarmonicToTime, PE->Type, Num_Element, PE->NbrNodes, case FORMAT_TIME_TABLE:
PE->x, PE->y, PE->z, Dummy, PE->Value) ; case FORMAT_SIMPLE_SPACE_TABLE:
break ; case FORMAT_VALUE_ONLY:
case FORMAT_SPACE_TABLE :
case FORMAT_TIME_TABLE :
case FORMAT_SIMPLE_SPACE_TABLE :
case FORMAT_VALUE_ONLY :
Tabular_PrintElement(PSO_P, PSO_P->Format, Time, TimeStep, NbTimeStep, Tabular_PrintElement(PSO_P, PSO_P->Format, Time, TimeStep, NbTimeStep,
NbrHarmonics, HarmonicToTime, PE->Type, Num_Element, NbrHarmonics, HarmonicToTime, PE->Type, Num_Element,
PE->NbrNodes, PE->x, PE->y, PE->z, Dummy, PE->Value) ; PE->NbrNodes, PE->x, PE->y, PE->z, Dummy, PE->Value);
break ; break;
case FORMAT_NODE_TABLE : case FORMAT_NODE_TABLE:
NodeTable_PrintElement(TimeStep, NbTimeStep, NbrHarmonics, PE); NodeTable_PrintElement(TimeStep, NbTimeStep, NbrHarmonics, PE);
break; break;
case FORMAT_ELEMENT_TABLE : case FORMAT_ELEMENT_TABLE:
ElementTable_PrintElement(TimeStep, NbTimeStep, NbrHarmonics, PE); ElementTable_PrintElement(TimeStep, NbTimeStep, NbrHarmonics, PE);
break; break;
case FORMAT_LOOP_ERROR : case FORMAT_LOOP_ERROR: StorePostOpResult(NbrHarmonics, PE); break;
StorePostOpResult(NbrHarmonics, PE);
break;
case FORMAT_ADAPT: case FORMAT_ADAPT:
if(PostStream){ if(PostStream) {
if(Dummy[4]) fprintf(PostStream, "%d\n", (int)Dummy[4]) ; if(Dummy[4]) fprintf(PostStream, "%d\n", (int)Dummy[4]);
fprintf(PostStream, "%d %g %g %g\n", fprintf(PostStream, "%d %g %g %g\n", (int)Dummy[0], Dummy[1], Dummy[2],
(int)Dummy[0], Dummy[1], Dummy[2], Dummy[3]); Dummy[3]);
} }
break ; break;
default : default: Message::Error("Unknown format in Format_PostElement");
Message::Error("Unknown format in Format_PostElement");
} }
if (PE->NbrNodes == 1 && if(PE->NbrNodes == 1 && PSO_P->Format != FORMAT_NODE_TABLE &&
PSO_P->Format != FORMAT_NODE_TABLE && PSO_P->Format != FORMAT_ELEMENT_TABLE) {
PSO_P->Format != FORMAT_ELEMENT_TABLE){ if(PSO_P->SendToServer && strcmp(PSO_P->SendToServer, "No")) {
if(PSO_P->SendToServer && strcmp(PSO_P->SendToServer, "No")){
std::vector<double> v; std::vector<double> v;
Export_Value(&PE->Value[0], v, PSO_P->SendToServerList); Export_Value(&PE->Value[0], v, PSO_P->SendToServerList);
Message::AddOnelabNumberChoice(PSO_P->SendToServer, v, PSO_P->Color, Message::AddOnelabNumberChoice(PSO_P->SendToServer, v, PSO_P->Color,
PSO_P->Units, PSO_P->Label, PSO_P->Visible, PSO_P->Units, PSO_P->Label, PSO_P->Visible,
PSO_P->Closed); PSO_P->Closed);
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* F o r m a t _ P o s t V a l u e */ /* F o r m a t _ P o s t V a l u e */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Format_PostValue(struct PostQuantity *PQ_P, void Format_PostValue(struct PostQuantity *PQ_P, struct PostSubOperation *PSO_P,
struct PostSubOperation *PSO_P,
int Format, char *Comma, int Group_FunctionType, int Format, char *Comma, int Group_FunctionType,
int iTime, double Time, int NbrTimeStep, int iTime, double Time, int NbrTimeStep, int iRegion,
int iRegion, int numRegion, int NbrRegion, int numRegion, int NbrRegion, int NbrHarmonics,
int NbrHarmonics, int HarmonicToTime, int FourierTransform, int HarmonicToTime, int FourierTransform,
int Flag_NoNewLine, int Flag_NoNewLine, struct Value *Value)
struct Value * Value) {
{ static int Size;
static int Size ; int j, k;
int j, k ; double TimeMH, Freq;
double TimeMH, Freq ; double x, y, z;
double x, y, z ; int flag_storeAllTimeResults, indexInTmpValues;
int flag_storeAllTimeResults, indexInTmpValues ; static struct Value TmpValue, *TmpValues;
static struct Value TmpValue, *TmpValues ; static double *Times;
static double *Times ;
flag_storeAllTimeResults = FourierTransform || PSO_P->TimeToHarmonic ; flag_storeAllTimeResults = FourierTransform || PSO_P->TimeToHarmonic;
indexInTmpValues = flag_storeAllTimeResults? iTime * NbrRegion : 0 ; indexInTmpValues = flag_storeAllTimeResults ? iTime * NbrRegion : 0;
if(iRegion == 0){ if(iRegion == 0) {
switch(Value->Type){ switch(Value->Type) {
case SCALAR : Size = 1 ; break ; case SCALAR: Size = 1; break;
case VECTOR : Size = 3 ; break ; case VECTOR: Size = 3; break;
case TENSOR_DIAG : Size = 3 ; break ; case TENSOR_DIAG: Size = 3; break;
case TENSOR_SYM : Size = 6 ; break ; case TENSOR_SYM: Size = 6; break;
case TENSOR : Size = 9 ; break ; case TENSOR: Size = 9; break;
} }
} }
if (Format == FORMAT_REGION_TABLE) { if(Format == FORMAT_REGION_TABLE) {
if(iRegion == 0){ if(iRegion == 0) {
if(PostStream == stdout || PostStream == stderr) { if(PostStream == stdout || PostStream == stderr) {
Message::Direct("%d%s", NbrRegion, Comma ? Comma : ""); Message::Direct("%d%s", NbrRegion, Comma ? Comma : "");
} }
else if(PostStream) { else if(PostStream) {
fprintf(PostStream, "%d%s\n", NbrRegion, Comma ? Comma : "") ; fprintf(PostStream, "%d%s\n", NbrRegion, Comma ? Comma : "");
} }
} }
std::ostringstream sstream; std::ostringstream sstream;
sstream.precision(16); sstream.precision(16);
sstream << numRegion; sstream << numRegion;
for (k = 0 ; k < NbrHarmonics ; k++) { for(k = 0; k < NbrHarmonics; k++) {
for(j = 0 ; j < Size ; j++) { for(j = 0; j < Size; j++) {
sstream << " " << Value->Val[MAX_DIM*k+j] ; sstream << " " << Value->Val[MAX_DIM * k + j];
if (Comma) sstream << Comma; if(Comma) sstream << Comma;
} }
} }
if(PostStream == stdout || PostStream == stderr) if(PostStream == stdout || PostStream == stderr)
Message::Direct(sstream.str().c_str()); Message::Direct(sstream.str().c_str());
else if(PostStream) else if(PostStream)
fprintf(PostStream, "%s\n", sstream.str().c_str()) ; fprintf(PostStream, "%s\n", sstream.str().c_str());
} }
else if (Format == FORMAT_GETDP) { else if(Format == FORMAT_GETDP) {
std::ostringstream sstream; std::ostringstream sstream;
sstream.precision(16); sstream.precision(16);
sstream << (PSO_P->ValueName? PSO_P->ValueName : PQ_P->Name); sstream << (PSO_P->ValueName ? PSO_P->ValueName : PQ_P->Name);
if (numRegion != NO_REGION) sstream << "~{" << numRegion << "}"; if(numRegion != NO_REGION) sstream << "~{" << numRegion << "}";
if (NbrHarmonics >= 2 || Size > 1) sstream << "() = {"; if(NbrHarmonics >= 2 || Size > 1)
else sstream << " ="; sstream << "() = {";
else
sstream << " =";
for (k = 0 ; k < NbrHarmonics ; k++) { for(k = 0; k < NbrHarmonics; k++) {
for(j = 0 ; j < Size ; j++) { for(j = 0; j < Size; j++) {
if (k || j) sstream << ","; if(k || j) sstream << ",";
sstream << " " << Value->Val[MAX_DIM*k+j] ; sstream << " " << Value->Val[MAX_DIM * k + j];
} }
} }
if (NbrHarmonics >= 2 || Size > 1) sstream << " };"; if(NbrHarmonics >= 2 || Size > 1) sstream << " };";
if(PostStream == stdout || PostStream == stderr) if(PostStream == stdout || PostStream == stderr)
Message::Direct(sstream.str().c_str()); Message::Direct(sstream.str().c_str());
else if(PostStream) else if(PostStream)
fprintf(PostStream, "%s\n", sstream.str().c_str()) ; fprintf(PostStream, "%s\n", sstream.str().c_str());
} }
else if (Format == FORMAT_GMSH && Flag_GMSH_VERSION != 2) { else if(Format == FORMAT_GMSH && Flag_GMSH_VERSION != 2) {
if (Group_FunctionType == NODESOF) if(Group_FunctionType == NODESOF)
Geo_GetNodesCoordinates(1, &numRegion, &x, &y, &z) ; Geo_GetNodesCoordinates(1, &numRegion, &x, &y, &z);
else { else {
x = y = z = 0.; x = y = z = 0.;
Message::Warning("Post Format \'Gmsh\' not adapted for global quantities s Message::Warning(
upported" "Post Format \'Gmsh\' not adapted for global quantities supported"
" by Regions. Zero coordinates are considered.") ; " by Regions. Zero coordinates are considered.");
} }
GmshParsed_PrintElement(Time, 0, 1, NbrHarmonics, HarmonicToTime, GmshParsed_PrintElement(Time, 0, 1, NbrHarmonics, HarmonicToTime,
POINT_ELEMENT, 1, &x, &y, &z, POINT_ELEMENT, 1, &x, &y, &z, Value);
Value) ;
} }
else if (Format == FORMAT_NXUNV) { else if(Format == FORMAT_NXUNV) {
if(PostStream) if(PostStream)
Unv_PrintRegion(PostStream, Comma ? 1 : 0, numRegion, NbrHarmonics, Unv_PrintRegion(PostStream, Comma ? 1 : 0, numRegion, NbrHarmonics, Size,
Size, Value, NXUnv_UnitFactor); Value, NXUnv_UnitFactor);
} }
else if (Format == FORMAT_LOOP_ERROR) { else if(Format == FORMAT_LOOP_ERROR) {
StorePostOpResult(NbrHarmonics, Value); StorePostOpResult(NbrHarmonics, Value);
} }
else if (Format == FORMAT_NODE_TABLE) { else if(Format == FORMAT_NODE_TABLE) {
// FIXME: this leads to output files without the total number of nodes at // FIXME: this leads to output files without the total number of nodes at
// the beginning (i.e. not compatible with NodeTable obtained e.g. for // the beginning (i.e. not compatible with NodeTable obtained e.g. for
// OnElementsOf) // OnElementsOf)
fprintf(PostStream, "%d", numRegion) ; fprintf(PostStream, "%d", numRegion);
Geo_GetNodesCoordinates(1, &numRegion, &x, &y, &z) ; Geo_GetNodesCoordinates(1, &numRegion, &x, &y, &z);
fprintf(PostStream, " %.16g %.16g %.16g", x,y,z) ; fprintf(PostStream, " %.16g %.16g %.16g", x, y, z);
for (k = 0 ; k < NbrHarmonics ; k++) { for(k = 0; k < NbrHarmonics; k++) {
for(j = 0 ; j < Size ; j++) { for(j = 0; j < Size; j++) {
fprintf(PostStream, " %.16g", Value->Val[MAX_DIM*k+j]) ; fprintf(PostStream, " %.16g", Value->Val[MAX_DIM * k + j]);
} }
} }
fprintf(PostStream, "\n") ; fprintf(PostStream, "\n");
} }
// else, for other FORMATs, e.g., FORMAT_FREQUENCY_TABLE // else, for other FORMATs, e.g., FORMAT_FREQUENCY_TABLE
else { else {
if(iRegion == 0){ if(iRegion == 0) {
if (!flag_storeAllTimeResults) if(!flag_storeAllTimeResults)
TmpValues = (struct Value*) Malloc(NbrRegion*sizeof(struct Value)) ; TmpValues = (struct Value *)Malloc(NbrRegion * sizeof(struct Value));
else{ else {
if (iTime == 0){ if(iTime == 0) {
TmpValues = (struct Value*) Malloc(NbrTimeStep*NbrRegion*sizeof(struct TmpValues = (struct Value *)Malloc(NbrTimeStep * NbrRegion *
Value)) ; sizeof(struct Value));
Times = (double*) Malloc(NbrTimeStep*sizeof(double)) ; Times = (double *)Malloc(NbrTimeStep * sizeof(double));
} }
Times[iTime] = Time ; Times[iTime] = Time;
} }
} }
Cal_CopyValue(Value, &TmpValues[indexInTmpValues+iRegion]) ; Cal_CopyValue(Value, &TmpValues[indexInTmpValues + iRegion]);
if (!flag_storeAllTimeResults && iRegion == NbrRegion-1) {
if (PostStream && HarmonicToTime == 1) { if(!flag_storeAllTimeResults && iRegion == NbrRegion - 1) {
switch (Format) { if(PostStream && HarmonicToTime == 1) {
case FORMAT_FREQUENCY_TABLE : switch(Format) {
case FORMAT_FREQUENCY_REGION_VALUE : case FORMAT_FREQUENCY_TABLE:
if (NbrHarmonics == 1){ case FORMAT_FREQUENCY_REGION_VALUE:
Message::Error("FrequencyTable format not allowed (only one harmonic) if(NbrHarmonics == 1) {
") ; Message::Error(
"FrequencyTable format not allowed (only one harmonic)");
return; return;
} }
break ;
case FORMAT_VALUE_ONLY :
break; break;
default : case FORMAT_VALUE_ONLY: break;
fprintf(PostStream, " %.16g", Time) ; default:
if (Comma) fprintf(PostStream, "%s", Comma); fprintf(PostStream, " %.16g", Time);
break ; if(Comma) fprintf(PostStream, "%s", Comma);
} break;
for (iRegion = 0 ; iRegion < NbrRegion ; iRegion++) { }
for (k = 0 ; k < NbrHarmonics ; k++) { for(iRegion = 0; iRegion < NbrRegion; iRegion++) {
if ((Format == FORMAT_FREQUENCY_TABLE || for(k = 0; k < NbrHarmonics; k++) {
Format == FORMAT_FREQUENCY_REGION_VALUE) if((Format == FORMAT_FREQUENCY_TABLE ||
&& !(k % 2) && iRegion == 0) { Format == FORMAT_FREQUENCY_REGION_VALUE) &&
Freq = Current.DofData->Val_Pulsation[0] / TWO_PI ; !(k % 2) && iRegion == 0) {
fprintf(PostStream, " %.16g", Freq) ; Freq = Current.DofData->Val_Pulsation[0] / TWO_PI;
if (Comma) fprintf(PostStream, "%s", Comma); fprintf(PostStream, " %.16g", Freq);
} if(Comma) fprintf(PostStream, "%s", Comma);
for(j = 0 ; j < Size ; j++) { }
if (Format != FORMAT_REGION_VALUE && for(j = 0; j < Size; j++) {
Format != FORMAT_FREQUENCY_REGION_VALUE) { if(Format != FORMAT_REGION_VALUE &&
fprintf(PostStream, " %.16g", Format != FORMAT_FREQUENCY_REGION_VALUE) {
TmpValues[indexInTmpValues+iRegion].Val[MAX_DIM*k+j]) ; fprintf(
if (Comma) fprintf(PostStream, "%s", Comma); PostStream, " %.16g",
TmpValues[indexInTmpValues + iRegion].Val[MAX_DIM * k + j]);
if(Comma) fprintf(PostStream, "%s", Comma);
} }
} }
} }
} }
if (Flag_NoNewLine || if(Flag_NoNewLine || Format == FORMAT_REGION_VALUE ||
Format == FORMAT_REGION_VALUE || Format == FORMAT_FREQUENCY_REGION_V Format == FORMAT_FREQUENCY_REGION_VALUE)
ALUE) fprintf(PostStream, " ");
fprintf(PostStream, " ") ; else
else fprintf(PostStream, "\n");
fprintf(PostStream, "\n") ; }
} else if(PostStream) {
else if(PostStream){ for(k = 0; k < HarmonicToTime; k++) {
for(k = 0 ; k < HarmonicToTime ; k++) { for(iRegion = 0; iRegion < NbrRegion; iRegion++) {
for (iRegion = 0 ; iRegion < NbrRegion ; iRegion++) { F_MHToTime0(k + iRegion, &TmpValues[indexInTmpValues + iRegion],
F_MHToTime0(k+iRegion, &TmpValues[indexInTmpValues+iRegion], &TmpValu &TmpValue, k, HarmonicToTime, &TimeMH);
e, if(iRegion == 0) {
k, HarmonicToTime, &TimeMH) ; fprintf(PostStream, " %.16g", TimeMH);
if (iRegion == 0) { if(Comma) fprintf(PostStream, "%s", Comma);
fprintf(PostStream, " %.16g", TimeMH) ;
if (Comma) fprintf(PostStream, "%s", Comma);
} }
for(j = 0 ; j < Size ; j++) { for(j = 0; j < Size; j++) {
fprintf(PostStream, " %.16g", TmpValue.Val[j]) ; fprintf(PostStream, " %.16g", TmpValue.Val[j]);
if (Comma) fprintf(PostStream, "%s", Comma); if(Comma) fprintf(PostStream, "%s", Comma);
} }
} }
fprintf(PostStream, "\n") ; fprintf(PostStream, "\n");
} }
} }
if (flag_storeAllTimeResults) Free(Times) ; if(flag_storeAllTimeResults) Free(Times);
Free(TmpValues) ; Free(TmpValues);
} }
else if (flag_storeAllTimeResults && else if(flag_storeAllTimeResults && iTime == NbrTimeStep - 1 &&
iTime == NbrTimeStep-1 && iRegion == NbrRegion-1) { iRegion == NbrRegion - 1) {
Pos_FourierTransform(NbrTimeStep, NbrRegion, Times, TmpValues, Size, 1,
Pos_FourierTransform(NbrTimeStep, NbrRegion, Times, TmpValues, Size, PSO_P->TimeToHarmonic, NULL, NULL, NULL);
1, PSO_P->TimeToHarmonic, NULL, NULL, NULL);
Free(Times); Free(Times);
Free(TmpValues) ; Free(TmpValues);
} }
} }
} }
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
/* P o s _ F o u r i e r T r a n s f o r m */ /* P o s _ F o u r i e r T r a n s f o r m */
/* ------------------------------------------------------------------------ */ /* ------------------------------------------------------------------------ */
void Pos_FourierTransform(int NbrTimeStep, int NbrRegion, void Pos_FourierTransform(int NbrTimeStep, int NbrRegion, double *Times,
double *Times, struct Value *TmpValues, int Size, struct Value *TmpValues, int Size, int TypeOutput,
int TypeOutput, int Nb_Freq_Select_0, int Nb_Freq_Select_0, int *NbrFreq,
int *NbrFreq, double **Frequencies, struct Value **Out double **Frequencies, struct Value **OutValues)
Values)
{ {
#if NEW_CODE #if NEW_CODE
*NbrFreq = (NbrTimeStep-1)/2+1; *NbrFreq = (NbrTimeStep - 1) / 2 + 1;
*Frequencies = (double *)Malloc(*NbrFreq*sizeof(double)); *Frequencies = (double *)Malloc(*NbrFreq * sizeof(double));
*OutValues = (struct Value *)Malloc(*NbrFreq*sizeof(struct Value)); *OutValues = (struct Value *)Malloc(*NbrFreq * sizeof(struct Value));
int nfft = *NbrFreq; int nfft = *NbrFreq;
kissfft<double> fft(nfft, false); kissfft<double> fft(nfft, false);
std::vector<std::complex<double> > inbuf(nfft); std::vector<std::complex<double> > inbuf(nfft);
std::vector<std::complex<double> > outbuf(nfft); std::vector<std::complex<double> > outbuf(nfft);
for (int k = 0; k < nfft; ++k) for(int k = 0; k < nfft; ++k)
inbuf[k]= std::complex<double>(rand()/(double)RAND_MAX - .5, inbuf[k] = std::complex<double>(rand() / (double)RAND_MAX - .5,
rand()/(double)RAND_MAX - .5); rand() / (double)RAND_MAX - .5);
fft.transform(&inbuf[0], &outbuf[0]); fft.transform(&inbuf[0], &outbuf[0]);
#else #else
int iTime, iRegion, k_fc, i_k, j, k; int iTime, iRegion, k_fc, i_k, j, k;
int N, Nhalf, NbrFourierComps; int N, Nhalf, NbrFourierComps;
double *val_FourierComps; double *val_FourierComps;
double val, val_r, val_i, norm, Period, w, v_cos, v_sin; double val, val_r, val_i, norm, Period, w, v_cos, v_sin;
N = NbrTimeStep-1; N = NbrTimeStep - 1;
Nhalf = N/2; Nhalf = N / 2;
// Nhalf = 2; // Nhalf = 2;
NbrFourierComps = Nhalf*2; NbrFourierComps = Nhalf * 2;
Period = Times[NbrTimeStep-1] - Times[0]; Period = Times[NbrTimeStep - 1] - Times[0];
w = TWO_PI/Period; w = TWO_PI / Period;
val_FourierComps = (double*) Malloc(NbrFourierComps*MAX_DIM*2*sizeof(double)) val_FourierComps =
; (double *)Malloc(NbrFourierComps * MAX_DIM * 2 * sizeof(double));
for (k_fc=-Nhalf; k_fc<Nhalf; k_fc++){ for(k_fc = -Nhalf; k_fc < Nhalf; k_fc++) {
i_k = Nhalf+k_fc; i_k = Nhalf + k_fc;
for (k=0; k<2; k++){ for(k = 0; k < 2; k++) {
for (j = 0 ; j < Size ; j++){ for(j = 0; j < Size; j++) {
val_FourierComps[(2*i_k+k)*MAX_DIM + j] = 0.; val_FourierComps[(2 * i_k + k) * MAX_DIM + j] = 0.;
} }
} }
} }
for (iTime=0; iTime<N; iTime++){ for(iTime = 0; iTime < N; iTime++) {
iRegion = 0; // only for 1 region now! iRegion = 0; // only for 1 region now!
for (k_fc=-Nhalf; k_fc<Nhalf; k_fc++){ for(k_fc = -Nhalf; k_fc < Nhalf; k_fc++) {
i_k = Nhalf+k_fc; i_k = Nhalf + k_fc;
v_cos = cos( k_fc*w*Times[iTime]); v_cos = cos(k_fc * w * Times[iTime]);
v_sin = sin(-k_fc*w*Times[iTime]); v_sin = sin(-k_fc * w * Times[iTime]);
for (j = 0 ; j < Size ; j++){ for(j = 0; j < Size; j++) {
val = TmpValues[iTime*NbrRegion+iRegion].Val[j]; val = TmpValues[iTime * NbrRegion + iRegion].Val[j];
val_FourierComps[(2*i_k+0)*MAX_DIM + j] += val * v_cos; val_FourierComps[(2 * i_k + 0) * MAX_DIM + j] += val * v_cos;
val_FourierComps[(2*i_k+1)*MAX_DIM + j] += val * v_sin; val_FourierComps[(2 * i_k + 1) * MAX_DIM + j] += val * v_sin;
} }
} }
} }
for (k_fc=-Nhalf; k_fc<Nhalf; k_fc++){ for(k_fc = -Nhalf; k_fc < Nhalf; k_fc++) {
i_k = Nhalf+k_fc; i_k = Nhalf + k_fc;
for (k=0; k<2; k++){ for(k = 0; k < 2; k++) {
for (j = 0 ; j < Size ; j++){ for(j = 0; j < Size; j++) {
val_FourierComps[(2*i_k+k)*MAX_DIM + j] /= N; val_FourierComps[(2 * i_k + k) * MAX_DIM + j] /= N;
} }
} }
} }
if(!PostStream) TypeOutput = 2; if(!PostStream) TypeOutput = 2;
int Nb_Freq_Select = (Nb_Freq_Select_0>0)? Nb_Freq_Select_0 : Nhalf-1; int Nb_Freq_Select = (Nb_Freq_Select_0 > 0) ? Nb_Freq_Select_0 : Nhalf - 1;
// Limited to cosine transform now // Limited to cosine transform now
if (TypeOutput == 2){ if(TypeOutput == 2) {
*NbrFreq = Nhalf+1; *NbrFreq = Nhalf + 1;
*Frequencies = (double *)Malloc(*NbrFreq*sizeof(double)); *Frequencies = (double *)Malloc(*NbrFreq * sizeof(double));
*OutValues = (struct Value *)Malloc(*NbrFreq*sizeof(struct Value)); *OutValues = (struct Value *)Malloc(*NbrFreq * sizeof(struct Value));
} }
// for (k_fc=-Nhalf; k_fc<Nhalf; k_fc++){ // for (k_fc=-Nhalf; k_fc<Nhalf; k_fc++){
for (k_fc=0; k_fc<=Nb_Freq_Select; k_fc++){ for(k_fc = 0; k_fc <= Nb_Freq_Select; k_fc++) {
i_k = Nhalf+k_fc; i_k = Nhalf + k_fc;
if (TypeOutput == 1) { if(TypeOutput == 1) { fprintf(PostStream, "%.16g", k_fc * w / TWO_PI); }
fprintf(PostStream, "%.16g", k_fc*w/TWO_PI) ; else if(TypeOutput == 2) {
} (*Frequencies)[k_fc] = k_fc * w / TWO_PI;
else if (TypeOutput == 2) {
(*Frequencies)[k_fc] = k_fc*w/TWO_PI;
(*OutValues)[k_fc].Type = TmpValues[0].Type; (*OutValues)[k_fc].Type = TmpValues[0].Type;
} }
for (k=0; k<2; k++){ for(k = 0; k < 2; k++) {
for(j = 0 ; j < Size ; j++){ for(j = 0; j < Size; j++) {
/* /*
if (k_fc != 0) if (k_fc != 0)
val = ((k==0)?1:-1) * 2 * val_FourierComps[(2*i_k+k)*MAX_DIM + j]; val = ((k==0)?1:-1) * 2 * val_FourierComps[(2*i_k+k)*MAX_DIM + j];
else else
val = val_FourierComps[(2*i_k+k)*MAX_DIM + j]; val = val_FourierComps[(2*i_k+k)*MAX_DIM + j];
*/ */
if (k_fc != 0) { if(k_fc != 0) {
val_r = 2 * val_FourierComps[(2*i_k+0)*MAX_DIM + j]; val_r = 2 * val_FourierComps[(2 * i_k + 0) * MAX_DIM + j];
val_i = -2 * val_FourierComps[(2*i_k+1)*MAX_DIM + j]; val_i = -2 * val_FourierComps[(2 * i_k + 1) * MAX_DIM + j];
norm = sqrt(SQU(val_r) + SQU(val_i)); norm = sqrt(SQU(val_r) + SQU(val_i));
// val = (k==0)? norm : -asin(val_i/norm); // Phase for CosineTransfor // val = (k==0)? norm : -asin(val_i/norm); // Phase for
m // CosineTransform
val = (k==0)? norm : atan2(val_i,val_r); // Phase for FourierTransform val =
(k == 0) ? norm : atan2(val_i, val_r); // Phase for FourierTransform
} }
else { else {
val = (k==0)? val_FourierComps[(2*i_k+k)*MAX_DIM + j] : 0.; val = (k == 0) ? val_FourierComps[(2 * i_k + k) * MAX_DIM + j] : 0.;
}
if (TypeOutput == 1) {
fprintf(PostStream, " %.16g", val);
}
if (k == 0 && TypeOutput == 2) {
(*OutValues)[k_fc].Val[j] = val;
} }
if(TypeOutput == 1) { fprintf(PostStream, " %.16g", val); }
if(k == 0 && TypeOutput == 2) { (*OutValues)[k_fc].Val[j] = val; }
} }
} }
if (TypeOutput == 1) { if(TypeOutput == 1) { fprintf(PostStream, "\n"); }
fprintf(PostStream, "\n") ;
}
} }
Free(val_FourierComps); Free(val_FourierComps);
#endif #endif
} }
 End of changes. 288 change blocks. 
1224 lines changed or deleted 1493 lines changed or added

Home  |  About  |  Features  |  All  |  Newest  |  Dox  |  Diffs  |  RSS Feeds  |  Screenshots  |  Comments  |  Imprint  |  Privacy  |  HTTP(S)