pos_format.c

cfd求解器使用与gmsh网格的求解

#define RCSID "$Id: Pos_Format.c,v 1.49 2006/02/26 00:42:59 geuzaine Exp $"
/*
 * Copyright (C) 1997-2006 P. Dular, C. Geuzaine
 *
 * This program is free software; you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation; either version 2 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
 * USA.
 *
 * Please report all bugs and problems to <[e-mail]>.
 */

#include "GetDP.h"
#include "GetDPVersion.h"
#include "GeoData.h"
#include "Data_Passive.h"
#include "Data_DefineE.h"
#include "Numeric.h"
#include "CurrentData.h"
#include "Magic.h"
#include "Pos_Iso.h"
#include "Pos_Format.h"
#include "Pos_Element.h"
#include "Pos_Formulation.h"
#include "F_Function.h"
#include "Cal_Value.h"
#include "Cal_Quantity.h"

#define NBR_MAX_ISO  200

/* bricolage: en attendant de trouver une meilleure solution pour les
   sorties ascii en format Gmsh (au niveau allocation memoire), je
   laisse la sortie GmshParsed par defaut. */

static List_T *PostElement_L = NULL ;
static List_T *TimeValue_L = NULL ;

static int     Gmsh_StartNewView = 0 ;
static int     NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT, NbTT;
static int     NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH;
static int     NbSI, NbVI, NbTI, NbSY, NbVY, NbTY;
static int     NbT2;
static List_T *SP = NULL, *VP = NULL, *TP = NULL;
static List_T *SL = NULL, *VL = NULL, *TL = NULL;
static List_T *ST = NULL, *VT = NULL, *TT = NULL;
static List_T *SQ = NULL, *VQ = NULL, *TQ = NULL;
static List_T *SS = NULL, *VS = NULL, *TS = NULL;
static List_T *SH = NULL, *VH = NULL, *TH = NULL;
static List_T *SI = NULL, *VI = NULL, *TI = NULL;
static List_T *SY = NULL, *VY = NULL, *TY = NULL;
static List_T *T2D = NULL, *T2C = NULL;

/* ------------------------------------------------------------------------ */
/*  F o r m a t _ P o s t F o r m a t / H e a d e r / F o o t e r           */
/* ------------------------------------------------------------------------ */

void  Format_PostFormat(int Format){
  
  GetDP_Begin("Format_PostFormat");

  switch(Format){
  case FORMAT_GMSH :
    if(Flag_BIN){/* bricolage */
      fprintf(PostStream, "$PostFormat /* Gmsh 1.2, %s */\n",
	      Flag_BIN ? "binary" : "ascii") ;
      fprintf(PostStream, "1.2 %d %d\n", Flag_BIN, (int)sizeof(double)) ;
      fprintf(PostStream, "$EndPostFormat\n") ;
    }
    break ;
  case FORMAT_GNUPLOT :
    fprintf(PostStream, "# GetDP %s, %s\n", GETDP_VERSION, 
	    Flag_BIN ? "binary" : "ascii") ;
    break ;
  }
  
  GetDP_End ;
}

void  Format_PostHeader(int Format, int Contour, 
			int NbTimeStep, int HarmonicToTime,
			int Type, int Order,
			char *Name1, char *Name2){

  char name[MAX_STRING_LENGTH] ;

  GetDP_Begin("Format_PostHeader");

  if(Contour){
    if(!PostElement_L) 
      PostElement_L = List_Create(20, 20, sizeof(struct PostElement*));
    else
      List_Reset(PostElement_L);
  }

  if(Name1 && Name2) {
    strcpy(name, Order ? Name1 : Name2) ;
    strcat(name, Get_StringForDefine(PostSubOperation_CombinationType, Type)) ;
    strcat(name, Order ? Name2 : Name1) ;
  }
  else if(Name1)
    strcpy(name, Name1) ;
  else if(Name2)
    strcpy(name, Name2) ;
  else
    strcpy(name, "unnamed");

  switch(Format){
  case FORMAT_GMSH_PARSED :
    fprintf(PostStream, "View \"%s\" {\n", name) ;
    Gmsh_StartNewView = 1 ;
    break ;
  case FORMAT_GMSH :
    if(Flag_BIN){ /* bricolage */
      fprintf(PostStream, "$View /* %s */\n", name);
      fprintf(PostStream, "%s ", name);
    }
    else{
      fprintf(PostStream, "View \"%s\" {\n", name) ;
    }
    Gmsh_StartNewView = 1 ;
    break ;
  case FORMAT_GNUPLOT :
    fprintf(PostStream, "# PostData '%s'\n", name);
    fprintf(PostStream, "# Type Num  X Y Z  N1 N2 N3  Values  <Values>...\n");
    break ;
  case FORMAT_ADAPT :
    fprintf(PostStream, "$Adapt /* %s */\n", name) ;
    break ;
  }

  GetDP_End ;
}

void Gmsh_StringStart(int Format, double x, double y, double style){
  double d;
  if(Flag_BIN){ /* bricolage: should use Format instead */
    List_Add(T2D, &x);
    List_Add(T2D, &y);
    List_Add(T2D, &style); 
    d = List_Nbr(T2C);
    List_Add(T2D, &d); 
    NbT2++;
  }
  else{
    fprintf(PostStream, "T2(%g,%g,%g){", x, y, style);
  }
}

void Gmsh_StringAdd(int Format, int first, char *text){
  int i;
  if(Flag_BIN){ /* bricolage: should use Format instead */
    for(i = 0; i < (int)strlen(text)+1; i++) 
      List_Add(T2C, &text[i]);
  }
  else{
    if(!first) 
      fprintf(PostStream, ",");
    fprintf(PostStream, "\"%s\"", text);
  }
}

void Gmsh_StringEnd(int Format){
  if(Flag_BIN){ /* bricolage: should use Format instead */
  }
  else{
    fprintf(PostStream, "};\n") ;
  }
}

void  Format_PostFooter(struct PostSubOperation *PSO_P, int Store){
  List_T  * Iso_L[NBR_MAX_ISO] ;
  double    IsoMin = 1.e200, IsoMax = -1.e200, IsoVal = 0.0, freq, valr, vali ;
  int       NbrIso = 0 ; 
  int       iPost, iNode, iIso, f, iTime, One=1, i, j, NbTimeStep ;
  char      tmp[1024];
  struct PostElement *PE ;
  
  GetDP_Begin("Format_PostFooter");
  
  if ( (PSO_P->Format == FORMAT_GMSH || PSO_P->Format == FORMAT_GMSH_PARSED) ){
    
    if( !(NbTimeStep = List_Nbr(PSO_P->TimeStep_L)) )
      NbTimeStep = List_Nbr(Current.DofData->Solutions);
    
    switch(PSO_P->Legend){
      
    case LEGEND_TIME:
      Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0],
		       PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]);
      for (i = 0 ; i < NbTimeStep ; i++) {
	Pos_InitAllSolutions(PSO_P->TimeStep_L, i) ;
	valr = Current.DofData->CurrentSolution->Time ;
	for (j = 0 ; j < Current.NbrHar ; j++){
	  sprintf(tmp, "Step %d/%d: Time = %g", i+1, NbTimeStep, valr);
	  Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp);
	}
      }
      Gmsh_StringEnd(PSO_P->Format);
      break;
      
    case LEGEND_FREQUENCY:
      if(Current.NbrHar > 1) {
	Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0],
			 PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]);
	for (i = 0 ; i < NbTimeStep ; i++) {
	  Pos_InitAllSolutions(PSO_P->TimeStep_L, i) ;
	  for (j = 0 ; j < Current.NbrHar ; j+=2) { 
	    freq = 0.5*Current.DofData->Val_Pulsation[j/2]/M_PI ;
	    sprintf(tmp, "%g Hz (Real Part: COSINUS)", freq);
	    Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp);
	    sprintf(tmp, "%g Hz (Imaginary Part: -SINUS)", freq);
	    Gmsh_StringAdd(PSO_P->Format, 0, tmp);
	  } 
	}
	Gmsh_StringEnd(PSO_P->Format);
      }
      break;
      
    case LEGEND_EIGENVALUES:
      Gmsh_StringStart(PSO_P->Format, PSO_P->LegendPosition[0],
		       PSO_P->LegendPosition[1], PSO_P->LegendPosition[2]);
      for (i = 0 ; i < NbTimeStep ; i++) {
	Pos_InitAllSolutions(PSO_P->TimeStep_L, i) ;
	valr = Current.DofData->CurrentSolution->Time ;
	vali = Current.DofData->CurrentSolution->TimeImag ;
	for (j = 0 ; j < Current.NbrHar ; j+=2) { 
	  sprintf(tmp, "Eigenvalue %d/%d: %g %s i * %g (Real Part)", 
		  i+1, NbTimeStep, valr, (vali > 0) ? "+" : "-", 
		  (vali > 0) ? vali : -vali);
	  Gmsh_StringAdd(PSO_P->Format, (!i && !j), tmp);
	  sprintf(tmp, "Eigenvalue %d/%d: %g %s i * %g (Imaginary Part)", 
		  i+1, NbTimeStep, valr, (vali > 0) ? "+" : "-", 
		  (vali > 0) ? vali : -vali);
	  Gmsh_StringAdd(PSO_P->Format, 0, tmp);
	}
      }
      Gmsh_StringEnd(PSO_P->Format);
      break;
    }
  }

  if(PSO_P->Iso){

    for(iPost = 0 ; iPost < List_Nbr(PostElement_L) ; iPost++){ 
      PE = *(struct PostElement**)List_Pointer(PostElement_L, iPost);
      for (iNode = 0 ; iNode < PE->NbrNodes ; iNode++ ){
	IsoMin = MIN(IsoMin, PE->Value[iNode].Val[0]) ;
	IsoMax = MAX(IsoMax, PE->Value[iNode].Val[0]) ;
      }
    }

    if((NbrIso = PSO_P->Iso) < 0)
      NbrIso = List_Nbr(PSO_P->Iso_L) ;
    if(NbrIso > NBR_MAX_ISO) 
      Msg(GERROR, "Too many Iso values");

    if(PSO_P->Format == FORMAT_GNUPLOT)
      fprintf(PostStream, "# NbIso = %d, Min = %g, Max = %g\n",
	      NbrIso, IsoMin, IsoMax) ;

    for(iIso = 0 ; iIso < NbrIso ; iIso++)
      Iso_L[iIso] = List_Create(10, 10, sizeof(struct PostElement*)) ;

    for(iPost = 0 ; iPost < List_Nbr(PostElement_L) ; iPost++){ 
      PE = *(struct PostElement**)List_Pointer(PostElement_L, iPost);
      for(iIso = 0 ; iIso < NbrIso ; iIso++){
	if(PSO_P->Iso > 0){
	  Cal_Iso(PE, Iso_L[iIso], IsoMin+iIso*(IsoMax-IsoMin)/(double)(NbrIso-1), 
		  IsoMin, IsoMax, PSO_P->DecomposeInSimplex) ;
	}
	else{
	  List_Read(PSO_P->Iso_L, iIso, &IsoVal) ;
	  Cal_Iso(PE, Iso_L[iIso], IsoVal, IsoMin, IsoMax, PSO_P->DecomposeInSimplex) ;
	}
      }
      if(!Store) Destroy_PostElement(PE);
    }

    for(iIso = 0 ; iIso < NbrIso ; iIso++){
      for(iPost = 0 ; iPost < List_Nbr(Iso_L[iIso]) ; iPost++){
	PE = *(struct PostElement**)List_Pointer(Iso_L[iIso], iPost) ;
	Format_PostElement(PSO_P->Format, 0, 0, 
			   Current.Time, 0, 1, 
			   Current.NbrHar, PSO_P->HarmonicToTime, 
			   NULL, PE,
			   PSO_P->ChangeOfCoordinates,
			   PSO_P->ChangeOfValues);
	Destroy_PostElement(PE) ;
      }
      List_Delete(Iso_L[iIso]) ;
      if(PSO_P->Format == FORMAT_GNUPLOT) fprintf(PostStream, "\n") ;
    }
  }

  switch(PSO_P->Format){
  case FORMAT_GMSH_PARSED :
    if(List_Nbr(TimeValue_L) > 1){
      fprintf(PostStream, "TIME{");
      for(iTime = 0; iTime < List_Nbr(TimeValue_L); iTime++){
	if(iTime) fprintf(PostStream, ",");
	fprintf(PostStream, "%g", *(double*)List_Pointer(TimeValue_L, iTime));
      }
      fprintf(PostStream, "};\n");
    }
    fprintf(PostStream, "};\n") ;
    break ;
  case FORMAT_GMSH :
    if(Flag_BIN){ /* bricolage */
      fprintf(PostStream, "%d %d %d %d %d %d %d %d %d %d %d %d %d %d %d %d "
	      "%d %d %d %d %d %d %d %d %d %d %d 0 0\n",
	      List_Nbr(TimeValue_L),
	      NbSP, NbVP, NbTP, NbSL, NbVL, NbTL, NbST, NbVT, NbTT, 
	      NbSQ, NbVQ, NbTQ, NbSS, NbVS, NbTS, NbSH, NbVH, NbTH,
	      NbSI, NbVI, NbTI, NbSY, NbVY, NbTY, NbT2, List_Nbr(T2C));
      if(Flag_BIN){
	f = LIST_FORMAT_BINARY;
	fwrite(&One, sizeof(int), 1, PostStream);
      }
      else{
	f = LIST_FORMAT_ASCII;
      }
      List_WriteToFile(TimeValue_L, PostStream, f); 
      List_WriteToFile(SP, PostStream, f); List_WriteToFile(VP, PostStream, f);
      List_WriteToFile(TP, PostStream, f); 
      List_WriteToFile(SL, PostStream, f); List_WriteToFile(VL, PostStream, f); 
      List_WriteToFile(TL, PostStream, f);
      List_WriteToFile(ST, PostStream, f); List_WriteToFile(VT, PostStream, f);
      List_WriteToFile(TT, PostStream, f); 
      List_WriteToFile(SQ, PostStream, f); List_WriteToFile(VQ, PostStream, f);
      List_WriteToFile(TQ, PostStream, f); 
      List_WriteToFile(SS, PostStream, f); List_WriteToFile(VS, PostStream, f);
      List_WriteToFile(TS, PostStream, f);
      List_WriteToFile(SH, PostStream, f); List_WriteToFile(VH, PostStream, f);
      List_WriteToFile(TH, PostStream, f);
      List_WriteToFile(SI, PostStream, f); List_WriteToFile(VI, PostStream, f);
      List_WriteToFile(TI, PostStream, f);
      List_WriteToFile(SY, PostStream, f); List_WriteToFile(VY, PostStream, f);
      List_WriteToFile(TY, PostStream, f);
      List_WriteToFile(T2D, PostStream, f); List_WriteToFile(T2C, PostStream, f);
      fprintf(PostStream, "\n");
      fprintf(PostStream, "$EndView\n");
    }
    else{
      if(List_Nbr(TimeValue_L) > 1){
	fprintf(PostStream, "TIME{");
	for(iTime = 0; iTime < List_Nbr(TimeValue_L); iTime++){
	  if(iTime) fprintf(PostStream, ",");
	  fprintf(PostStream, "%g", *(double*)List_Pointer(TimeValue_L, iTime));
	}
	fprintf(PostStream, "};\n");
      }
      fprintf(PostStream, "};\n") ;
    }
    break ;
  case FORMAT_ADAPT :
    fprintf(PostStream, "$EndAdapt\n");
    break ;
  }

  GetDP_End ;
}

/* ------------------------------------------------------------------------ */
/*  F o r m a t _ G m s h P a r s e d                                       */
/* ------------------------------------------------------------------------ */

void  Format_GmshParsed(double Time, int TimeStep, int NbTimeStep, int NbHarmonic, 
			int HarmonicToTime, int Type, int NbrNodes, 
			double *x, double *y, double *z, struct Value *Value){
  int     i,j,k;