[cig-commits] r13676 - in long/3D/Gale/trunk: . src/PICellerator/Utils/src

walter at geodynamics.org walter at geodynamics.org
Sat Dec 13 02:16:26 PST 2008


Author: walter
Date: 2008-12-13 02:16:26 -0800 (Sat, 13 Dec 2008)
New Revision: 13676

Added:
   long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.c
   long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.h
   long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.meta
Modified:
   long/3D/Gale/trunk/
   long/3D/Gale/trunk/src/PICellerator/Utils/src/Init.c
   long/3D/Gale/trunk/src/PICellerator/Utils/src/SConscript
   long/3D/Gale/trunk/src/PICellerator/Utils/src/Utils.h
   long/3D/Gale/trunk/src/PICellerator/Utils/src/types.h
Log:
 r2429 at dante:  boo | 2008-12-13 01:50:57 -0800
 Add HydrostaticTerm for the analytic pressure and density functions



Property changes on: long/3D/Gale/trunk
___________________________________________________________________
Name: svk:merge
   - 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2428
   + 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2429

Added: long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.c
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.c	                        (rev 0)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.c	2008-12-13 10:16:26 UTC (rev 13676)
@@ -0,0 +1,375 @@
+/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+**
+** Copyright (C), 2003-2006, Victorian Partnership for Advanced Computing (VPAC) Ltd, 110 Victoria Street,
+**	Melbourne, 3053, Australia.
+** Copyright (c) 2005-2006, Monash Cluster Computing, Building 28, Monash University Clayton Campus,
+**	Victoria, 3800, Australia
+**
+** Primary Contributing Organisations:
+**	Victorian Partnership for Advanced Computing Ltd, Computational Software Development - http://csd.vpac.org
+**	Australian Computational Earth Systems Simulator - http://www.access.edu.au
+**	Monash Cluster Computing - http://www.mcc.monash.edu.au
+**
+** Contributors:
+**	Robert Turnbull, Research Assistant, Monash University. (robert.turnbull at sci.monash.edu.au)
+**	Patrick D. Sunter, Software Engineer, VPAC. (patrick at vpac.org)
+**	Alan H. Lo, Computational Engineer, VPAC. (alan at vpac.org)
+**	Stevan M. Quenette, Senior Software Engineer, VPAC. (steve at vpac.org)
+**	David May, PhD Student, Monash University (david.may at sci.monash.edu.au)
+**	Vincent Lemiale, Postdoctoral Fellow, Monash University. (vincent.lemiale at sci.monash.edu.au)
+**	Julian Giordani, Research Assistant, Monash University. (julian.giordani at sci.monash.edu.au)
+**	Louis Moresi, Associate Professor, Monash University. (louis.moresi at sci.monash.edu.au)
+**	Luke J. Hodkinson, Computational Engineer, VPAC. (lhodkins at vpac.org)
+**	Raquibul Hassan, Computational Engineer, VPAC. (raq at vpac.org)
+**	David Stegman, Postdoctoral Fellow, Monash University. (david.stegman at sci.monash.edu.au)
+**	Wendy Sharples, PhD Student, Monash University (wendy.sharples at sci.monash.edu.au)
+**
+**  This library is free software; you can redistribute it and/or
+**  modify it under the terms of the GNU Lesser General Public
+**  License as published by the Free Software Foundation; either
+**  version 2.1 of the License, or (at your option) any later version.
+**
+**  This library 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
+**  Lesser General Public License for more details.
+**
+**  You should have received a copy of the GNU Lesser General Public
+**  License along with this library; if not, write to the Free Software
+**  Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
+**
+** $Id: /cig/src/PICellerator/Utils/src/HydrostaticTerm.c 2420 2008-12-12T20:12:30.951338Z boo  $
+**
+**~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
+
+
+#include <mpi.h>
+#include <StGermain/StGermain.h>
+#include <StgFEM/StgFEM.h>
+#include <PICellerator/Voronoi/Voronoi.h>
+#include <PICellerator/PopulationControl/PopulationControl.h>
+#include <PICellerator/Weights/Weights.h>
+#include <PICellerator/MaterialPoints/MaterialPoints.h>
+
+#include "types.h"
+#include "HydrostaticTerm.h"
+#include "MaterialSwarmVariable.h"
+
+#include <assert.h>
+#include <string.h>
+
+/* Textual name of this class */
+const Type HydrostaticTerm_Type = "HydrostaticTerm";
+
+HydrostaticTerm* HydrostaticTerm_New(Name name,
+                                     double upper_density,
+                                     double upper_alpha,
+                                     double lower_density,
+                                     double lower_alpha,
+                                     double height,
+                                     double material_boundary,
+                                     double T_0,
+                                     double linear_coefficient,
+                                     double exponential_coefficient1,
+                                     double exponential_coefficient2,
+                                     double gravity)
+{
+	HydrostaticTerm* self = (HydrostaticTerm*) _HydrostaticTerm_DefaultNew( name );
+
+	HydrostaticTerm_InitAll(self,
+                                upper_density,
+                                upper_alpha,
+                                lower_density,
+                                lower_alpha,
+                                height,
+                                material_boundary,
+                                T_0,
+                                linear_coefficient,
+                                exponential_coefficient1,
+                                exponential_coefficient2,
+                                gravity);
+	return self;
+}
+
+/* Creation implementation / Virtual constructor */
+HydrostaticTerm* _HydrostaticTerm_New( 
+		SizeT                                               sizeOfSelf,  
+		Type                                                type,
+		Stg_Class_DeleteFunction*                           _delete,
+		Stg_Class_PrintFunction*                            _print,
+		Stg_Class_CopyFunction*                             _copy, 
+		Stg_Component_DefaultConstructorFunction*           _defaultConstructor,
+		Stg_Component_ConstructFunction*                    _construct,
+		Stg_Component_BuildFunction*                        _build,
+		Stg_Component_InitialiseFunction*                   _initialise,
+		Stg_Component_ExecuteFunction*                      _execute,
+		Stg_Component_DestroyFunction*                      _destroy,
+                double upper_density,
+                double upper_alpha,
+                double lower_density,
+                double lower_alpha,
+                double height,
+                double material_boundary,
+                double T_0,
+                double linear_coefficient,
+                double exponential_coefficient1,
+                double exponential_coefficient2,
+                double gravity,
+		Name                                                name )
+{
+  HydrostaticTerm* self;
+  
+  /* Allocate memory */
+  assert( sizeOfSelf >= sizeof(HydrostaticTerm) );
+  self = (ForceTerm*)_Stg_Component_New(sizeOfSelf,
+                                        type, 
+                                        _delete,
+                                        _print,
+                                        _copy,
+                                        _defaultConstructor,
+                                        _construct,
+                                        _build,
+                                        _initialise,
+                                        _execute,
+                                        _destroy,
+                                        name,
+                                        NON_GLOBAL );
+
+  /* Virtual info */
+  
+  self->upper_density=upper_density;
+  self->upper_alpha=upper_alpha;
+  self->lower_density=lower_density;
+  self->lower_alpha=lower_alpha;
+  self->height=height;
+  self->material_boundary=material_boundary;
+  self->T_0=T_0;
+  self->linear_coefficient=linear_coefficient;
+  self->exponential_coefficient1=exponential_coefficient1;
+  self->exponential_coefficient2=exponential_coefficient2;
+  self->gravity=gravity;
+  return self;
+}
+
+void _HydrostaticTerm_Init(HydrostaticTerm* self,
+                           double upper_density,
+                           double upper_alpha,
+                           double lower_density,
+                           double lower_alpha,
+                           double height,
+                           double material_boundary,
+                           double T_0,
+                           double linear_coefficient,
+                           double exponential_coefficient1,
+                           double exponential_coefficient2,
+                           double gravity)
+{
+  self->isConstructed    = True;
+
+  self->upper_density=upper_density;
+  self->upper_alpha=upper_alpha;
+  self->lower_density=lower_density;
+  self->lower_alpha=lower_alpha;
+  self->height=height;
+  self->material_boundary=material_boundary;
+  self->T_0=T_0;
+  self->linear_coefficient=linear_coefficient;
+  self->exponential_coefficient1=exponential_coefficient1;
+  self->exponential_coefficient2=exponential_coefficient2;
+  self->gravity=gravity;
+}
+
+void HydrostaticTerm_InitAll(void* forceTerm,
+                             double upper_density,
+                             double upper_alpha,
+                             double lower_density,
+                             double lower_alpha,
+                             double height,
+                             double material_boundary,
+                             double T_0,
+                             double linear_coefficient,
+                             double exponential_coefficient1,
+                             double exponential_coefficient2,
+                             double gravity)
+{
+  HydrostaticTerm* self = (HydrostaticTerm*) forceTerm;
+  _HydrostaticTerm_Init( self,
+                         upper_density,
+                         upper_alpha,
+                         lower_density,
+                         lower_alpha,
+                         height,
+                         material_boundary,
+                         T_0,
+                         linear_coefficient,
+                         exponential_coefficient1,
+                         exponential_coefficient2,
+                         gravity);
+}
+
+void _HydrostaticTerm_Delete( void* forceTerm ) {
+}
+
+void _HydrostaticTerm_Print( void* forceTerm, Stream* stream ) {
+}
+
+void* _HydrostaticTerm_DefaultNew( Name name ) {
+	return (void*)_HydrostaticTerm_New( 
+		sizeof(HydrostaticTerm), 
+		HydrostaticTerm_Type,
+		_HydrostaticTerm_Delete,
+		_HydrostaticTerm_Print,
+		NULL,
+		_HydrostaticTerm_DefaultNew,
+		_HydrostaticTerm_Construct,
+		_HydrostaticTerm_Build,
+		_HydrostaticTerm_Initialise,
+		_HydrostaticTerm_Execute,
+		_HydrostaticTerm_Destroy,
+                0,0,0,0,0,0,0,0,0,0,0,
+		name );
+}
+
+void _HydrostaticTerm_Construct( void* forceTerm, Stg_ComponentFactory* cf,
+                                 void* data ) {
+  HydrostaticTerm* self = (HydrostaticTerm*)forceTerm;
+  double upper_density,upper_alpha,lower_density,lower_alpha,height,
+    material_boundary,T_0,linear_coefficient,exponential_coefficient1,
+    exponential_coefficient2,gravity;
+
+  upper_density=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"upperDensity",0.0);
+  upper_alpha=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"upperAlpha",0.0);
+  lower_density=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"lowerDensity",0.0);
+  lower_alpha=Stg_ComponentFactory_GetDouble(cf,self->name,"lowerAlpha",0.0);
+  height=Stg_ComponentFactory_GetDouble(cf,self->name,"height",0.0);
+  material_boundary=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"materialBoundary",0.0);
+  T_0=Stg_ComponentFactory_GetDouble(cf,self->name,"T_0",0.0);
+  linear_coefficient=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"linearCoefficient",0.0);
+  exponential_coefficient1=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"exponentialCoefficient1",0.0);
+  exponential_coefficient2=
+    Stg_ComponentFactory_GetDouble(cf,self->name,"exponentialCoefficient2",0.0);
+  gravity=Stg_ComponentFactory_GetDouble(cf,self->name,"gravity",0.0);
+
+  _HydrostaticTerm_Init( self, upper_density,
+                         upper_alpha,
+                         lower_density,
+                         lower_alpha,
+                         height,
+                         material_boundary,
+                         T_0,
+                         linear_coefficient,
+                         exponential_coefficient1,
+                         exponential_coefficient2,
+                         gravity );
+}
+
+void _HydrostaticTerm_Build( void* forceTerm, void* data ) {
+}
+
+void _HydrostaticTerm_Initialise( void* forceTerm, void* data ) {
+}
+
+void _HydrostaticTerm_Execute( void* forceTerm, void* data ) {
+}
+
+void _HydrostaticTerm_Destroy( void* forceTerm, void* data ) {
+}
+
+double HydrostaticTerm_Density( void* forceTerm, Coord coord)
+{
+  HydrostaticTerm *self=(HydrostaticTerm *)forceTerm;
+  double h=self->height-coord[1];
+  double T, density, alpha;
+
+  /* First, get the temperature */
+  if(h<0)
+    {
+      T=self->T_0;
+    }
+  else
+    {
+      T=self->T_0 + self->linear_coefficient*h
+        + self->exponential_coefficient1
+        *(1-exp(-self->exponential_coefficient2*h));
+    }
+
+  /* Then, get the density */
+  if(coord[1]>self->height)
+    {
+      density=0;
+    }
+  else if(coord[1]>self->material_boundary)
+    {
+      density=self->upper_density;;
+      alpha=self->upper_alpha;
+    }
+  else
+    {
+      density=self->lower_density;
+      alpha=self->lower_alpha;
+    }
+
+  return density*(1-alpha*T);
+}
+
+double HydrostaticTerm_Pressure_Analytic(double density, double gravity,
+                                         double h, double alpha,
+                                         double T_0, double A,
+                                         double B, double C)
+{
+  return density*gravity*h*(1-alpha*(T_0 + A*h/2 + B*(1+exp(-C*h)/(C*h))));
+}
+
+double HydrostaticTerm_Pressure( void* forceTerm, Coord coord)
+{
+  HydrostaticTerm *self=(HydrostaticTerm *)forceTerm;
+  double h=self->height-coord[1];
+  double p;
+
+  if(coord[1]>self->height)
+    {
+      p=0;
+    }
+  else if(coord[1]>self->material_boundary)
+    {
+      p=HydrostaticTerm_Pressure_Analytic(self->upper_density, self->gravity, h,
+                                          self->upper_alpha,
+                                          self->T_0,
+                                          self->linear_coefficient,
+                                          self->exponential_coefficient1,
+                                          self->exponential_coefficient2);
+    }
+  else
+    {
+      p=HydrostaticTerm_Pressure_Analytic(self->upper_density,
+                                          self->gravity,
+                                          self->height-self->material_boundary,
+                                          self->upper_alpha,
+                                          self->T_0,
+                                          self->linear_coefficient,
+                                          self->exponential_coefficient1,
+                                          self->exponential_coefficient2);
+      - HydrostaticTerm_Pressure_Analytic(self->lower_density, self->gravity,
+                                          self->height-self->material_boundary,
+                                          self->lower_alpha,
+                                          self->T_0,
+                                          self->linear_coefficient,
+                                          self->exponential_coefficient1,
+                                          self->exponential_coefficient2);
+        + HydrostaticTerm_Pressure_Analytic(self->lower_density, self->gravity,
+                                            h,
+                                            self->lower_alpha,
+                                            self->T_0,
+                                            self->linear_coefficient,
+                                            self->exponential_coefficient1,
+                                            self->exponential_coefficient2);
+    }
+  return p;
+}
+

Added: long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.h
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.h	                        (rev 0)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.h	2008-12-13 10:16:26 UTC (rev 13676)
@@ -0,0 +1,137 @@
+/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+**
+** Copyright (C), 2003-2006, Victorian Partnership for Advanced Computing (VPAC) Ltd, 110 Victoria Street,
+**	Melbourne, 3053, Australia.
+** Copyright (c) 2005-2006, Monash Cluster Computing, Building 28, Monash University Clayton Campus,
+**	Victoria, 3800, Australia
+**
+** Primary Contributing Organisations:
+**	Victorian Partnership for Advanced Computing Ltd, Computational Software Development - http://csd.vpac.org
+**	Australian Computational Earth Systems Simulator - http://www.access.edu.au
+**	Monash Cluster Computing - http://www.mcc.monash.edu.au
+**
+** Contributors:
+**	Robert Turnbull, Research Assistant, Monash University. (robert.turnbull at sci.monash.edu.au)
+**	Patrick D. Sunter, Software Engineer, VPAC. (patrick at vpac.org)
+**	Alan H. Lo, Computational Engineer, VPAC. (alan at vpac.org)
+**	Stevan M. Quenette, Senior Software Engineer, VPAC. (steve at vpac.org)
+**	David May, PhD Student, Monash University (david.may at sci.monash.edu.au)
+**	Vincent Lemiale, Postdoctoral Fellow, Monash University. (vincent.lemiale at sci.monash.edu.au)
+**	Julian Giordani, Research Assistant, Monash University. (julian.giordani at sci.monash.edu.au)
+**	Louis Moresi, Associate Professor, Monash University. (louis.moresi at sci.monash.edu.au)
+**	Luke J. Hodkinson, Computational Engineer, VPAC. (lhodkins at vpac.org)
+**	Raquibul Hassan, Computational Engineer, VPAC. (raq at vpac.org)
+**	David Stegman, Postdoctoral Fellow, Monash University. (david.stegman at sci.monash.edu.au)
+**	Wendy Sharples, PhD Student, Monash University (wendy.sharples at sci.monash.edu.au)
+**
+**  This library is free software; you can redistribute it and/or
+**  modify it under the terms of the GNU Lesser General Public
+**  License as published by the Free Software Foundation; either
+**  version 2.1 of the License, or (at your option) any later version.
+**
+**  This library 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
+**  Lesser General Public License for more details.
+**
+**  You should have received a copy of the GNU Lesser General Public
+**  License along with this library; if not, write to the Free Software
+**  Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA  02110-1301  USA
+**
+*/
+
+
+#ifndef __PICellerator_Utils_HydrostaticTerm_h__
+#define __PICellerator_Utils_HydrostaticTerm_h__
+
+	/** Textual name of this class */
+	extern const Type HydrostaticTerm_Type;
+
+	/** HydrostaticTerm class contents */
+	#define __HydrostaticTerm \
+		/* General info */ \
+		__Stg_Component \
+		\
+		/* Virtual info */ \
+                double upper_density; \
+                double upper_alpha; \
+                double lower_density; \
+                double lower_alpha; \
+                double height; \
+                double material_boundary; \
+                double T_0; \
+                double linear_coefficient; \
+                double exponential_coefficient1; \
+                double exponential_coefficient2; \
+                double gravity; \
+
+
+	struct HydrostaticTerm { __HydrostaticTerm };
+
+	HydrostaticTerm* HydrostaticTerm_New( 
+		Name                                                name,
+                double upper_density,
+                double upper_alpha,
+                double lower_density,
+                double lower_alpha,
+                double height,
+                double material_boundary,
+                double T_0,
+                double linear_coefficient,
+                double exponential_coefficient1,
+                double exponential_coefficient2,
+                double gravity);
+
+	HydrostaticTerm* _HydrostaticTerm_New( 
+		SizeT                                               sizeOfSelf,  
+		Type                                                type,
+		Stg_Class_DeleteFunction*                           _delete,
+		Stg_Class_PrintFunction*                            _print,
+		Stg_Class_CopyFunction*                             _copy, 
+		Stg_Component_DefaultConstructorFunction*           _defaultConstructor,
+		Stg_Component_ConstructFunction*                    _construct,
+		Stg_Component_BuildFunction*                        _build,
+		Stg_Component_InitialiseFunction*                   _initialise,
+		Stg_Component_ExecuteFunction*                      _execute,
+		Stg_Component_DestroyFunction*                      _destroy,
+                double upper_density,
+                double upper_alpha,
+                double lower_density,
+                double lower_alpha,
+                double height,
+                double material_boundary,
+                double T_0,
+                double linear_coefficient,
+                double exponential_coefficient1,
+                double exponential_coefficient2,
+                double gravity,
+		Name                                                name );
+	
+        void HydrostaticTerm_InitAll(void* forceTerm,
+                                     double upper_density,
+                                     double upper_alpha,
+                                     double lower_density,
+                                     double lower_alpha,
+                                     double height,
+                                     double material_boundary,
+                                     double T_0,
+                                     double linear_coefficient,
+                                     double exponential_coefficient1,
+                                     double exponential_coefficient2,
+                                     double gravity);
+
+	void _HydrostaticTerm_Delete( void* forceTerm );
+	void _HydrostaticTerm_Print( void* forceTerm, Stream* stream );
+
+	void* _HydrostaticTerm_DefaultNew( Name name ) ;
+        void _HydrostaticTerm_Construct( void* forceTerm,
+                                         Stg_ComponentFactory* cf,
+                                         void* data ) ;
+	void _HydrostaticTerm_Build( void* forceTerm, void* data ) ;
+	void _HydrostaticTerm_Initialise( void* forceTerm, void* data ) ;
+	void _HydrostaticTerm_Execute( void* forceTerm, void* data ) ;
+	void _HydrostaticTerm_Destroy( void* forceTerm, void* data ) ;
+        double HydrostaticTerm_Density( void* forceTerm, Coord coord);
+        double HydrostaticTerm_Pressure( void* forceTerm, Coord coord);
+
+#endif

Added: long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.meta
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.meta	                        (rev 0)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/HydrostaticTerm.meta	2008-12-13 10:16:26 UTC (rev 13676)
@@ -0,0 +1,32 @@
+<?xml version="1.0"?>
+<!DOCTYPE StGermainData SYSTEM "stgermain.dtd">
+<StGermainData xmlns="http://www.vpac.org/StGermain/XML_IO_Handler/Jun2003">
+
+<param name="Name">HydrostaticTerm</param>
+<param name="Author">...</param>
+<param name="Organisation">VPAC and MCC</param>
+<param name="Project">PICellerator</param>
+<param name="Location">./PICellerator/Utils/src/</param>
+<param name="Project Web">https://csd.vpac.org/twiki/bin/view/PICellerator/WebHome</param>
+<param name="Copyright">Copyright (C) 2005 VPAC and Monash Cluster Computing.</param>
+<param name="License">https://csd.vpac.org/twiki/bin/view/Stgermain/SoftwareLicense http://www.opensource.org/licenses/bsd-license.php</param>
+<param name="Parent">ForceTerm</param>
+<param name="Reference">...</param>
+<param name="Summary">...</param>
+<param name="Description"></param>
+<param name="Equation"></param>
+
+<!--Now the interesting stuff-->
+
+
+<list name="Params">
+
+</list>
+
+<list name="Dependencies">
+</list>
+<!-- Add an exmaple XML if possible -->
+<param name="Example">
+</param>
+
+</StGermainData>

Modified: long/3D/Gale/trunk/src/PICellerator/Utils/src/Init.c
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/Init.c	2008-12-13 10:16:21 UTC (rev 13675)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/Init.c	2008-12-13 10:16:26 UTC (rev 13676)
@@ -62,9 +62,12 @@
 	Stg_ComponentRegister_Add( componentsRegister, BuoyancyForceTerm_Type,            "0", _BuoyancyForceTerm_DefaultNew );
 	Stg_ComponentRegister_Add( componentsRegister, BuoyancyForceTermThermoChem_Type,            "0", _BuoyancyForceTermThermoChem_DefaultNew );
 
+	Stg_ComponentRegister_Add( componentsRegister, HydrostaticTerm_Type,            "0", _HydrostaticTerm_DefaultNew );
+
 	RegisterParent( BuoyancyForceTerm_Type,     ForceTerm_Type );
 	RegisterParent( BuoyancyForceTermThermoChem_Type,     ForceTerm_Type );
 	RegisterParent( MaterialSwarmVariable_Type, SwarmVariable_Type );
 	RegisterParent( PCDVC_Type,                 DVCWeights_Type );
+	RegisterParent( HydrostaticTerm_Type,       Stg_Component_Type );
 	return True;
 }

Modified: long/3D/Gale/trunk/src/PICellerator/Utils/src/SConscript
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/SConscript	2008-12-13 10:16:21 UTC (rev 13675)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/SConscript	2008-12-13 10:16:26 UTC (rev 13676)
@@ -22,6 +22,7 @@
 header_groups['PICellerator/Utils']=Split("""BuoyancyForceTerm.h
 BuoyancyForceTermThermoChem.h
 Finalise.h
+HydrostaticTerm.h
 Init.h
 MaterialSwarmVariable.h
 PCDVC.h
@@ -32,6 +33,7 @@
 c_files=Split("""BuoyancyForceTerm.c
 BuoyancyForceTermThermoChem.c
 Finalise.c
+HydrostaticTerm.c
 Init.c
 PCDVC.c
 MaterialSwarmVariable.c""")
@@ -39,6 +41,7 @@
 
 meta_files=Split("""BuoyancyForceTerm.meta
 BuoyancyForceTermThermoChem.meta
+HydrostaticTerm.meta
 PCDVC.meta
 MaterialSwarmVariable.meta""")
 

Modified: long/3D/Gale/trunk/src/PICellerator/Utils/src/Utils.h
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/Utils.h	2008-12-13 10:16:21 UTC (rev 13675)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/Utils.h	2008-12-13 10:16:26 UTC (rev 13676)
@@ -48,6 +48,7 @@
 	#include "BuoyancyForceTerm.h"
 	#include "BuoyancyForceTermThermoChem.h"
 	#include "MaterialSwarmVariable.h"
+	#include "HydrostaticTerm.h"
 
 	#include "Init.h"
 	#include "Finalise.h"

Modified: long/3D/Gale/trunk/src/PICellerator/Utils/src/types.h
===================================================================
--- long/3D/Gale/trunk/src/PICellerator/Utils/src/types.h	2008-12-13 10:16:21 UTC (rev 13675)
+++ long/3D/Gale/trunk/src/PICellerator/Utils/src/types.h	2008-12-13 10:16:26 UTC (rev 13676)
@@ -57,6 +57,7 @@
 	typedef struct BuoyancyForceTerm                BuoyancyForceTerm;
 	typedef struct MaterialSwarmVariable            MaterialSwarmVariable;
 	typedef struct BuoyancyForceTermThermoChem      BuoyancyForceTermThermoChem;
+        typedef struct HydrostaticTerm      HydrostaticTerm;
         typedef struct PCDVC                        PCDVC;
 
 #endif 



More information about the CIG-COMMITS mailing list