pos_quantity.c

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

#define RCSID "$Id: Pos_Quantity.c,v 1.20 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 "Treatment_Formulation.h"
#include "Pos_Formulation.h"
#include "Pos_Quantity.h"
#include "Get_DofOfElement.h"
#include "GeoData.h"
#include "Cal_Quantity.h"
#include "Get_Geometry.h"
#include "Data_Passive.h"
#include "Data_DefineE.h"
#include "CurrentData.h"
#include "Tools.h"


/* ------------------------------------------------------------------------ */
/*  C a l _ P o s t Q u a n t i t y                                         */
/* ------------------------------------------------------------------------ */

void Cal_PostQuantity(struct PostQuantity    *PostQuantity_P, 
		      struct DefineQuantity  *DefineQuantity_P0,
		      struct QuantityStorage *QuantityStorage_P0,
		      List_T *Support_L,
		      struct Element         *Element, 
		      double u, double v, double w, 
		      struct Value *Value) {

  struct PostQuantityTerm  PostQuantityTerm ;

  List_T   *InRegion_L ;
  int       i_PQT, Type_Quantity ;

  GetDP_Begin("Cal_PostQuantity");

  /* mettre tout a zero: on ne connait pas a priori le type de retour */
  /* (default type and value returned if Type_Quantity == -1) */

  Cal_ZeroValue(Value);
  Value->Type = SCALAR; 

  /* Loop on PostQuantity Terms */
  /* ... with sum of results if common supports (In ...) */

  for (i_PQT = 0 ; i_PQT < List_Nbr(PostQuantity_P->PostQuantityTerm) ; i_PQT++) {
    
    List_Read(PostQuantity_P->PostQuantityTerm, i_PQT, &PostQuantityTerm) ;
    
    InRegion_L = (PostQuantityTerm.InIndex < 0)?  NULL :
      ((struct Group *)List_Pointer(Problem_S.Group, 
				    PostQuantityTerm.InIndex))->InitialList ;

    if (!Support_L)  Type_Quantity = PostQuantityTerm.Type ;
    else             Type_Quantity = GLOBALQUANTITY ; /* Always if Support */

    if (InRegion_L) {
      if (Element->Num != NO_ELEMENT) {
	if (!List_Search(InRegion_L, &Element->Region, fcmp_int)) { 
	  Type_Quantity = -1 ; 
	}
      }
      else {
	if (Type_Quantity == GLOBALQUANTITY) {
	  /* Plus de test ici... vu que le OnRegion de la PostOperation n'a rien
	     a voir avec le support d'une integration ...
	  if (!List_Search(InRegion_L, &Current.Region, fcmp_int)) {
	    Type_Quantity = -1 ;
	  }
	  */
	  /* Il faut plutot voir si il existe au moins une region de InRegion_L
	     qui soit dans Support_L ... cela est fait apres, pour chaque element */
	}
	else if (Type_Quantity != INTEGRALQUANTITY) {
	  Type_Quantity = -1 ;
	}
      }
    }
    /* else if !InRegion_L -> No filter, i.e. globally defined quantity */

    /* ---------------------------- */
    /* Local or Integral quantities */
    /* ---------------------------- */

    if (Type_Quantity == LOCALQUANTITY || Type_Quantity == INTEGRALQUANTITY) {

      Pos_LocalOrIntegralQuantity(PostQuantity_P,
				  DefineQuantity_P0, QuantityStorage_P0,
				  &PostQuantityTerm, Element, Type_Quantity,
				  u, v, w, Value) ;
    }

    /* ----------------- */
    /* Global quantities */
    /* ----------------- */

    else if (Type_Quantity == GLOBALQUANTITY) {

      Pos_GlobalQuantity(PostQuantity_P,
			 DefineQuantity_P0, QuantityStorage_P0,
			 &PostQuantityTerm, Element, InRegion_L, Support_L, Value) ;
    }
   
  }  /* for i_PQT ... */

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  P o s _ G l o b a l Q u a n t i t y                                     */
/* ------------------------------------------------------------------------ */

void Pos_GlobalQuantity(struct PostQuantity    *PostQuantity_P,
			struct DefineQuantity  *DefineQuantity_P0,
			struct QuantityStorage *QuantityStorage_P0,
			struct PostQuantityTerm  *PostQuantityTerm_P,
			struct Element         *ElementEmpty, 
			List_T  *InRegion_L, List_T  *Support_L,
			struct Value *Value) {

  struct DefineQuantity    *DefineQuantity_P ;
  struct QuantityStorage   *QuantityStorage_P ;
  struct FunctionSpace     *FunctionSpace_P ;
  struct GlobalQuantity    *GlobalQuantity_P ;
  struct Value              TermValue ;

  int  k, Index_DefineQuantity ;

  int    Nbr_Element, i_Element ;
  struct Element  Element ;
  int    Type_Quantity ;

  GetDP_Begin("Pos_GlobalQuantity");

  if (PostQuantityTerm_P->EvaluationType == LOCAL &&
      List_Search(InRegion_L, &Current.Region, fcmp_int)) {

    for (k = 0 ; k < PostQuantityTerm_P->NbrQuantityIndex ; k++) {	  
      Index_DefineQuantity = PostQuantityTerm_P->QuantityIndexTable[k] ;
      DefineQuantity_P     = DefineQuantity_P0  + Index_DefineQuantity ;
      QuantityStorage_P    = QuantityStorage_P0 + Index_DefineQuantity ;

      if (QuantityStorage_P->NumLastElementForFunctionSpace != Current.Region) {
	QuantityStorage_P->NumLastElementForFunctionSpace = Current.Region ;
	QuantityStorage_P->FunctionSpace = FunctionSpace_P =
	  (struct FunctionSpace*)
	  List_Pointer(Problem_S.FunctionSpace,
		       DefineQuantity_P->FunctionSpaceIndex) ;
	GlobalQuantity_P = (struct GlobalQuantity*)
	  List_Pointer
	  (QuantityStorage_P->FunctionSpace->GlobalQuantity,
	   *(int *)List_Pointer(DefineQuantity_P->IndexInFunctionSpace, 0)) ;

	if (DefineQuantity_P->Type == GLOBALQUANTITY)
	  Get_DofOfRegion(Current.Region, GlobalQuantity_P,
			  FunctionSpace_P, QuantityStorage_P) ;
      }
    }

    Cal_WholeQuantity
      (Current.Element = ElementEmpty,
       QuantityStorage_P0, PostQuantityTerm_P->WholeQuantity,
       Current.u = 0., Current.v = 0., Current.w = 0., -1, -1, &TermValue) ;

    Value->Type = TermValue.Type;
    Cal_AddValue(Value,&TermValue,Value);

  }  /* if LOCAL && ... */

  else if (PostQuantityTerm_P->EvaluationType == INTEGRAL) {

    Nbr_Element = Geo_GetNbrGeoElements() ;
    Get_InitDofOfElement(&Element) ;

    Type_Quantity = LOCALQUANTITY ; /* Attention... il faut se comprendre: */
    /* il s'agit de grandeurs locales qui seront integrees */
    for (i_Element = 0 ; i_Element < Nbr_Element; i_Element++) {
      Progress(i_Element, Nbr_Element, "Accumulate: ") ;

      Element.GeoElement = Geo_GetGeoElement(i_Element) ;
      Element.Num    = Element.GeoElement->Num ;
      Element.Type   = Element.GeoElement->Type ;
      Current.Region = Element.Region = Element.GeoElement->Region ;

      /* Filter: only elements in both InRegion_L and Support_L are considered */
      if ((!InRegion_L || List_Search(InRegion_L, &Element.Region, fcmp_int)) &&
	  (!Support_L || List_Search(Support_L, &Element.Region, fcmp_int))) {
	Get_NodesCoordinatesOfElement(&Element) ;
	Current.x = Element.x[0];
	Current.y = Element.y[0];
	Current.z = Element.z[0]; 
	Pos_LocalOrIntegralQuantity(PostQuantity_P,
				    DefineQuantity_P0, QuantityStorage_P0,
				    PostQuantityTerm_P, &Element, Type_Quantity,
				    0., 0., 0., Value) ;
      }
    }  /* for i_Element ... */
  }  /* if INTEGRAL ... */

  GetDP_End ;
}


/* ------------------------------------------------------------------------ */
/*  C a l _ P o s t C u m u l a t i v e Q u a n t i t y                     */
/* ------------------------------------------------------------------------ */

void Cal_PostCumulativeQuantity(List_T                 *Region_L,
				int                    SupportIndex,
				List_T                 *TimeStep_L, 
				struct PostQuantity    *PostQuantity_P, 
				struct DefineQuantity  *DefineQuantity_P0,
				struct QuantityStorage *QuantityStorage_P0,
				struct Value          **Values) {
  struct Element Element ;
  List_T *Support_L ;
  int i, NbrTimeStep ;

  GetDP_Begin("Cal_PostCumulativeQuantity");

  Support_L = ((struct Group *)
	       List_Pointer(Problem_S.Group, SupportIndex))->InitialList ;

  NbrTimeStep = List_Nbr(TimeStep_L) ;
  *Values = (struct Value *)Malloc(NbrTimeStep*sizeof(struct Value)) ;

  Element.Num = NO_ELEMENT ;

  for(i = 0 ; i < NbrTimeStep ; i++) {
    Pos_InitAllSolutions(TimeStep_L, i) ;