dofdata.c

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

#define RCSID "$Id: DofData.c,v 1.48 2006/03/13 20:49:51 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]>.
 *
 * Contributor(s):
 *   Johan Gyselinck
 *   Ruth Sabariego
 */

#include "GetDP.h"
#include "GetDPVersion.h"
#include "DofData.h"
#include "Tools.h"
#include "Magic.h"
#include "CurrentData.h"
#include "Numeric.h"

FILE  * File_PRE, * File_RES, * File_TMP ;

struct DofData  * CurrentDofData ;

int fcmp_Dof(const void * a, const void * b) {
  int Result ;

  if ((Result = ((struct Dof *)a)->NumType  - ((struct Dof *)b)->NumType) != 0)
    return Result ;
  if ((Result = ((struct Dof *)a)->Entity   - ((struct Dof *)b)->Entity)  != 0)
    return Result ;
  return        ((struct Dof *)a)->Harmonic - ((struct Dof *)b)->Harmonic ;
}

/* ------------------------------------------------------------------------ */
/*  D o f _ I n i t D o f D a t a                                           */
/* ------------------------------------------------------------------------ */

void  Dof_InitDofData(struct DofData * DofData_P, int Num,
		      int ResolutionIndex, int SystemIndex,
		      char * Name_SolverDataFile) {

  int  Index ;

  GetDP_Begin("Dof_InitDofData");

  DofData_P->Num = Num ;

  DofData_P->ResolutionIndex = ResolutionIndex ;
  DofData_P->SystemIndex     = SystemIndex ;

  DofData_P->FunctionSpaceIndex = NULL ;
  DofData_P->TimeFunctionIndex  = List_Create(10, 5, sizeof(int)) ;
  Index = 0 ;  List_Add(DofData_P->TimeFunctionIndex, &Index) ;
  DofData_P->Pulsation = NULL ;  DofData_P->Val_Pulsation = NULL ;
  DofData_P->NbrHar = 1 ;

  DofData_P->NbrAnyDof = 0 ;  DofData_P->NbrDof = 0 ;
  DofData_P->DofTree = Tree_Create(sizeof(struct Dof), fcmp_Dof) ;
  DofData_P->DofList = NULL ;

  DofData_P->NbrPart = 0 ;
  DofData_P->Nnz = NULL ;

  DofData_P->SolverDataFileName = Name_SolverDataFile ;
  DofData_P->Flag_Init[0] = 0 ;
  DofData_P->Flag_Init[1] = 0 ;
  DofData_P->Flag_Init[2] = 0 ;
  DofData_P->Flag_Init[3] = 0 ;

  DofData_P->Flag_Only = 0 ;
  DofData_P->Flag_InitOnly[0] = 0 ;
  DofData_P->Flag_InitOnly[1] = 0 ;
  DofData_P->Flag_InitOnly[2] = 0 ;

  DofData_P->OnlyTheseMatrices = NULL;

  DofData_P->Solutions = NULL ;
  DofData_P->CurrentSolution = NULL ;

  DofData_P->DummyDof = NULL ;

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  D o f _ S e t C u r r e n t D o f D a t a                               */
/* ------------------------------------------------------------------------ */

void  Dof_SetCurrentDofData(struct DofData * DofData_P) {

  GetDP_Begin("Dof_SetCurrentDofData");

  CurrentDofData = DofData_P ;

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  F i l e s   . . .                                                       */
/* ------------------------------------------------------------------------ */

/* ------------------------------------------------------------------------ */
/*  D o f _ O p e n F i l e                                                 */
/* ------------------------------------------------------------------------ */

void  Dof_OpenFile(int Type, char * Name, char * Mode) {
  char  * Extension, FileName[MAX_FILE_NAME_LENGTH] ;
  FILE  * File_X ;

  GetDP_Begin("Dof_OpenFile");

  switch (Type) {
  case DOF_PRE :  Extension = ".pre" ;  break ;
  case DOF_RES :  Extension = ""     ;  break ;
  case DOF_TMP :  Extension = ""     ;  break ;
  default      :  Extension = ".pre" ;  break ;
  }

  strcpy(FileName, Name) ; strcat(FileName, Extension) ;

  if (!(File_X = fopen(FileName, Mode)))
    Msg(GERROR,"Unable to open file '%s'", FileName) ;

  switch (Type) {
  case DOF_PRE :  File_PRE = File_X ;  break ;
  case DOF_RES :  File_RES = File_X ;  break ;
  case DOF_TMP :  File_TMP = File_RES ; File_RES = File_X ;  break ;
  default      :  break ;
  }

  GetDP_End ;
}

/* ------------------------------------------------------------------------ */
/*  D o f _ C l o s e F i l e                                               */
/* ------------------------------------------------------------------------ */

void  Dof_CloseFile(int Type) {

  GetDP_Begin("Dof_CloseFile");

  switch (Type) {
  case DOF_PRE :  fclose(File_PRE) ;  break ;
  case DOF_RES :  fclose(File_RES) ;  break ;
  case DOF_TMP :  fclose(File_RES) ;  File_RES = File_TMP ; break ;
  }

  GetDP_End ;
}

/* ------------------------------------------------------------------------ */
/*  D o f _ F l u s h F i l e                                               */
/* ------------------------------------------------------------------------ */

void  Dof_FlushFile(int Type) {

  GetDP_Begin("Dof_FlushFile");

  switch (Type) {
  case DOF_PRE :  fflush(File_PRE) ;  break ;
  case DOF_RES :  fflush(File_RES) ;  break ;
  }

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  D o f _ W r i t e F i l e P R E 0                                       */
/* ------------------------------------------------------------------------ */

void  Dof_WriteFilePRE0(int Num_Resolution, char * Name_Resolution,
			int Nbr_DofData) {

  GetDP_Begin("Dof_WriteFilePRE0");

  fprintf(File_PRE, "$Resolution /* '%s' */\n", Name_Resolution) ;
  fprintf(File_PRE, "%d %d\n", Num_Resolution, Nbr_DofData) ;
  fprintf(File_PRE, "$EndResolution\n") ;

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  D o f _ R e a d F i l e P R E 0                                         */
/* ------------------------------------------------------------------------ */

void  Dof_ReadFilePRE0(int * Num_Resolution, int * Nbr_DofData) {
  char  String[MAX_STRING_LENGTH] ;

  GetDP_Begin("Dof_ReadFilePRE0");

  do { 
    fgets(String, MAX_STRING_LENGTH, File_PRE) ; 
    if (feof(File_PRE))  break ;
  } while (String[0] != '$') ;

  if (feof(File_PRE)) Msg(GERROR, "$Resolution field not found in file");

  if (!strncmp(&String[1], "Resolution", 10)) {
    fscanf(File_PRE, "%d %d", Num_Resolution, Nbr_DofData) ;
  }

  do {
    fgets(String, MAX_STRING_LENGTH, File_PRE) ;
    if (feof(File_PRE)) Msg(GERROR, "Prematured end of file");
  } while (String[0] != '$') ;

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  D o f _ W r i t e F i l e P R E                                         */
/* ------------------------------------------------------------------------ */

static int * Nnz ;

void  Dof_WriteFilePRE(struct DofData * DofData_P) {

  int  i, Nbr_Index ;
  struct Dof  * Dof_P0 ;

  GetDP_Begin("Dof_WriteFilePRE");

  fprintf(File_PRE, "$DofData /* #%d */\n", DofData_P->Num) ;

  fprintf(File_PRE, "%d %d\n",
	  DofData_P->ResolutionIndex, DofData_P->SystemIndex) ;

  Nbr_Index = List_Nbr(DofData_P->FunctionSpaceIndex) ;
  fprintf(File_PRE, "%d", Nbr_Index) ;
  for (i = 0 ; i < Nbr_Index ; i++)
    fprintf(File_PRE, " %d",
	    *((int *)List_Pointer(DofData_P->FunctionSpaceIndex, i))) ;
  fprintf(File_PRE, "\n") ;

  Nbr_Index = List_Nbr(DofData_P->TimeFunctionIndex) ;
  fprintf(File_PRE, "%d", Nbr_Index) ;
  for (i = 0 ; i < Nbr_Index ; i++)
    fprintf(File_PRE, " %d",
	    *((int *)List_Pointer(DofData_P->TimeFunctionIndex, i))) ;
  fprintf(File_PRE, "\n") ;

  fprintf(File_PRE, "%d", DofData_P->NbrPart) ;
  for(i = 0 ; i < DofData_P->NbrPart+1 ; i++)
    fprintf(File_PRE, " %d", DofData_P->Part[i]) ;
  fprintf(File_PRE, "\n") ;

  fprintf(File_PRE, "%d %d\n",
	  (DofData_P->DofTree)? Tree_Nbr(DofData_P->DofTree) : DofData_P->NbrAnyDof,
	  DofData_P->NbrDof) ;

  Nnz = DofData_P->Nnz ;

  if (DofData_P->DofTree)
    Tree_Action(DofData_P->DofTree, Dof_WriteDofPRE) ;
  else {
    if (DofData_P->NbrAnyDof){
      Dof_P0 = (struct Dof *)List_Pointer(DofData_P->DofList, 0) ;
      for (i = 0 ; i < DofData_P->NbrAnyDof ; i++)
	Dof_WriteDofPRE(Dof_P0 + i, NULL) ;
    }
  }
  fprintf(File_PRE, "$EndDofData\n") ;
  fflush(File_PRE) ;

  GetDP_End ;
}

/* ------------------------------- */
/*  D o f _ W r i t e D o f P R E  */
/* ------------------------------- */

void  Dof_WriteDofPRE(void * a, void * b) {

  struct Dof  * Dof_P ;

  GetDP_Begin("Dof_WriteDofPRE");

  Dof_P = (struct Dof *) a ;

  fprintf(File_PRE, "%d %d %d %d ",
	  Dof_P->NumType, Dof_P->Entity, Dof_P->Harmonic, Dof_P->Type) ;

  switch (Dof_P->Type) {
  case DOF_UNKNOWN :
    fprintf(File_PRE, "%d %d\n", Dof_P->Case.Unknown.NumDof,
	    Nnz[Dof_P->Case.Unknown.NumDof-1]) ;
    break ;
  case DOF_FIXEDWITHASSOCIATE :
    fprintf(File_PRE, "%d ", Dof_P->Case.FixedAssociate.NumDof) ;
    LinAlg_PrintScalar(File_PRE, &Dof_P->Val);
    fprintf(File_PRE, " %d\n", Dof_P->Case.FixedAssociate.TimeFunctionIndex) ;
    break ;
  case DOF_FIXED :
    LinAlg_PrintScalar(File_PRE, &Dof_P->Val);
    fprintf(File_PRE, " %d\n", Dof_P->Case.FixedAssociate.TimeFunctionIndex) ;
    break ;
  case DOF_FIXED_SOLVE :
    fprintf(File_PRE, "%d\n", Dof_P->Case.FixedAssociate.TimeFunctionIndex) ;
    break ;
  case DOF_UNKNOWN_INIT :
    fprintf(File_PRE, "%d ", Dof_P->Case.Unknown.NumDof) ;
    LinAlg_PrintScalar(File_PRE, &Dof_P->Val);
    fprintf(File_PRE, " %d\n", Nnz[Dof_P->Case.Unknown.NumDof-1]) ;
    break ;
  case DOF_LINK :
    fprintf(File_PRE, "%.16g %d\n",
	    Dof_P->Case.Link.Coef, Dof_P->Case.Link.EntityRef) ;
    break ;
  case DOF_LINKCPLX :
    fprintf(File_PRE, "%.16g %.16g %d\n",
	    Dof_P->Case.Link.Coef, Dof_P->Case.Link.Coef2, 
	    Dof_P->Case.Link.EntityRef) ;
    break ;
  }

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  D o f _ R e a d F i l e P R E                                           */
/* ------------------------------------------------------------------------ */

void  Dof_ReadFilePRE(struct DofData * DofData_P) {

  int         i, Nbr_Index, Int ;
  struct Dof  Dof ;
  char        String[MAX_STRING_LENGTH] ;

  GetDP_Begin("Dof_ReadFilePRE");

  do { 
    fgets(String, MAX_STRING_LENGTH, File_PRE) ; 
    if (feof(File_PRE))  break ;
  } while (String[0] != '$') ;

  if (feof(File_PRE)) Msg(GERROR, "$DofData field not found in file");

  if (!strncmp(&String[1], "DofData", 7)) {

    fscanf(File_PRE, "%d %d",
	   &DofData_P->ResolutionIndex, &DofData_P->SystemIndex) ;

    fscanf(File_PRE, "%d", &Nbr_Index) ;
    DofData_P->FunctionSpaceIndex = List_Create(Nbr_Index, 1, sizeof(int)) ;
    for (i = 0 ; i < Nbr_Index ; i++) {
      fscanf(File_PRE, "%d", &Int) ;
      List_Add(DofData_P->FunctionSpaceIndex, &Int) ;
    }

    fscanf(File_PRE, "%d", &Nbr_Index) ;
    DofData_P->TimeFunctionIndex = List_Create(Nbr_Index, 1, sizeof(int)) ;
    for (i = 0 ; i < Nbr_Index ; i++) {
      fscanf(File_PRE, "%d", &Int) ;
      List_Add(DofData_P->TimeFunctionIndex, &Int) ;
    }

    fscanf(File_PRE, "%d", &DofData_P->NbrPart) ;
    for(i = 0 ; i < DofData_P->NbrPart+1 ; i++)
      fscanf(File_PRE, "%d", &DofData_P->Part[i]) ;

    fscanf(File_PRE, "%d %d", &DofData_P->NbrAnyDof, &DofData_P->NbrDof) ;

    DofData_P->DofList = List_Create(DofData_P->NbrAnyDof, 1, sizeof(struct Dof)) ;

    if(!DofData_P->Nnz)
      DofData_P->Nnz = (int*)Malloc(DofData_P->NbrDof * sizeof(int)) ;

    for (i = 0 ; i < DofData_P->NbrAnyDof ; i++) {

      fscanf(File_PRE, "%d %d %d %d",
	     &Dof.NumType, &Dof.Entity, &Dof.Harmonic, &Dof.Type) ;

      switch (Dof.Type) {
      case DOF_UNKNOWN :
	fscanf(File_PRE, "%d", &Dof.Case.Unknown.NumDof) ;
	fscanf(File_PRE, "%d", &DofData_P->Nnz[Dof.Case.Unknown.NumDof-1]) ;
	break ;
      case DOF_FIXEDWITHASSOCIATE :
	fscanf(File_PRE, "%d", &Dof.Case.FixedAssociate.NumDof) ;
	LinAlg_ScanScalar(File_PRE, &Dof.Val) ;
	fscanf(File_PRE, "%d", &Dof.Case.FixedAssociate.TimeFunctionIndex) ;
	break ;
      case DOF_FIXED :
	LinAlg_ScanScalar(File_PRE, &Dof.Val) ;
	fscanf(File_PRE, "%d", &Dof.Case.FixedAssociate.TimeFunctionIndex) ;
	break ;
      case DOF_FIXED_SOLVE :
	fscanf(File_PRE, "%d", &Dof.Case.FixedAssociate.TimeFunctionIndex) ;
	break ;
      case DOF_UNKNOWN_INIT :
	fscanf(File_PRE, "%d", &Dof.Case.Unknown.NumDof) ;
	LinAlg_ScanScalar(File_PRE, &Dof.Val) ;
	fscanf(File_PRE, "%d", &DofData_P->Nnz[Dof.Case.Unknown.NumDof-1]) ;
	break ;
      case DOF_LINK :
	fscanf(File_PRE, "%lf %d",
	       &Dof.Case.Link.Coef, &Dof.Case.Link.EntityRef) ;
	Dof.Case.Link.Dof = NULL ;
	break ;
      case DOF_LINKCPLX :
	fscanf(File_PRE, "%lf %lf %d",
	       &Dof.Case.Link.Coef, &Dof.Case.Link.Coef2, &Dof.Case.Link.EntityRef) ;
	Dof.Case.Link.Dof = NULL ;
	break ;
      }

      List_Add(DofData_P->DofList, &Dof) ;
    }
  }

  do {
    fgets(String, MAX_STRING_LENGTH, File_PRE) ;
    if (feof(File_PRE)) Msg(GERROR, "Prematured end of file");
  } while (String[0] != '$') ;

  Dof_InitDofType(DofData_P) ;

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  D o f _ W r i t e F i l e R E S 0                                       */
/* ------------------------------------------------------------------------ */

void  Dof_WriteFileRES0(char * Name_File, int Format) {

  GetDP_Begin("Dof_WriteFileRES0");

  LinAlg_SequentialBegin();

  Dof_OpenFile(DOF_RES, Name_File, (char*)(Format ? "wb" : "w")) ;
  fprintf(File_RES, "$ResFormat /* GetDP %s, %s */\n", GETDP_VERSION, 
	  Format ? "binary" : "ascii") ;
  fprintf(File_RES, "1.1 %d\n", Format) ;
  fprintf(File_RES, "$EndResFormat\n") ;