[cig-commits] r4787 - in long/3D/Gale/trunk: . src/Underworld/Rheology/src

walter at geodynamics.org walter at geodynamics.org
Wed Oct 11 11:39:24 PDT 2006


Author: walter
Date: 2006-10-11 11:39:23 -0700 (Wed, 11 Oct 2006)
New Revision: 4787

Added:
   long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.c
   long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.h
   long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.meta
Modified:
   long/3D/Gale/trunk/
   long/3D/Gale/trunk/src/Underworld/Rheology/src/Init.c
   long/3D/Gale/trunk/src/Underworld/Rheology/src/Rheology.h
   long/3D/Gale/trunk/src/Underworld/Rheology/src/SConscript
   long/3D/Gale/trunk/src/Underworld/Rheology/src/types.h
Log:
 r706 at earth:  boo | 2006-10-11 11:38:40 -0700
 Add MohrCoulomb



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

Modified: long/3D/Gale/trunk/src/Underworld/Rheology/src/Init.c
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/Init.c	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/Init.c	2006-10-11 18:39:23 UTC (rev 4787)
@@ -73,6 +73,7 @@
 	Stg_ComponentRegister_Add( componentRegister, Byerlee_Type,                 "0", _Byerlee_DefaultNew );
 	Stg_ComponentRegister_Add( componentRegister, DruckerPrager_Type,           "0", _DruckerPrager_DefaultNew );
 	Stg_ComponentRegister_Add( componentRegister, FaultingMoresiMuhlhaus2006_Type,             "0", _FaultingMoresiMuhlhaus2006_DefaultNew );
+	Stg_ComponentRegister_Add( componentRegister, MohrCoulomb_Type,             "0", _MohrCoulomb_DefaultNew );
 
 	Stg_ComponentRegister_Add( componentRegister, StrainWeakening_Type,         "0", _StrainWeakening_DefaultNew );
 	Stg_ComponentRegister_Add( componentRegister, BuiterStrainWeakening_Type,   "0", _BuiterStrainWeakening_DefaultNew );
@@ -102,6 +103,7 @@
 	RegisterParent( Byerlee_Type,                    VonMises_Type );
 	RegisterParent( DruckerPrager_Type,              VonMises_Type );
 	RegisterParent( FaultingMoresiMuhlhaus2006_Type, YieldRheology_Type );
+	RegisterParent( MohrCoulomb_Type, YieldRheology_Type );
 	
 	RegisterParent( StrainWeakening_Type,         TimeIntegratee_Type );
 	RegisterParent( BuiterStrainWeakening_Type,   StrainWeakening_Type );

Added: long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.c
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.c	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.c	2006-10-11 18:39:23 UTC (rev 4787)
@@ -0,0 +1,410 @@
+/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+** Copyright (c) 2005, Monash Cluster Computing 
+** All rights reserved.
+** Redistribution and use in source and binary forms, with or without modification,
+** are permitted provided that the following conditions are met:
+**
+** 		* Redistributions of source code must retain the above copyright notice, 
+** 			this list of conditions and the following disclaimer.
+** 		* Redistributions in binary form must reproduce the above copyright 
+**			notice, this list of conditions and the following disclaimer in the 
+**			documentation and/or other materials provided with the distribution.
+** 		* Neither the name of the Monash University nor the names of its contributors 
+**			may be used to endorse or promote products derived from this software 
+**			without specific prior written permission.
+**
+** THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 
+** "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, 
+** THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR 
+** PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS 
+** BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 
+** CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 
+** SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 
+** HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 
+** LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT 
+** OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+**
+**
+** Contact:
+*%		Louis Moresi - Louis.Moresi at sci.monash.edu.au
+*%
+** Contributors:
+*+		Robert Turnbull
+*+		Vincent Lemiale
+*+		Louis Moresi
+*+		David May
+*+		David Stegman
+*+		Mirko Velic
+*+		Patrick Sunter
+*+		Julian Giordani
+*+
+** $Id$
+** 
+**~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
+
+
+#include <mpi.h>
+#include <StGermain/StGermain.h>
+#include <StG_FEM/StG_FEM.h>
+#include <PICellerator/PICellerator.h>
+
+#include "types.h"
+#include "RheologyClass.h"
+#include "StrainWeakening.h"
+#include "YieldRheology.h"
+#include "MohrCoulomb.h"
+#include "ConstitutiveMatrix.h"
+
+#include <assert.h>
+#include <string.h>
+
+/* Textual name of this class - This is a global pointer which is used for times when you need to refer to class and not a particular instance of a class */
+const Type MohrCoulomb_Type = "MohrCoulomb";
+
+/* Private Constructor: This will accept all the virtual functions for this class as arguments. */
+MohrCoulomb* _MohrCoulomb_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,
+		Rheology_ModifyConstitutiveMatrixFunction*         _modifyConstitutiveMatrix,
+		YieldRheology_GetYieldCriterionFunction*           _getYieldCriterion,
+		YieldRheology_GetYieldIndicatorFunction*           _getYieldIndicator,
+		YieldRheology_HasYieldedFunction*                  _hasYielded,
+		Name                                               name ) 
+{
+	MohrCoulomb*					self;
+
+	/* Call private constructor of parent - this will set virtual functions of parent and continue up the hierarchy tree. At the beginning of the tree it will allocate memory of the size of object and initialise all the memory to zero. */
+	assert( sizeOfSelf >= sizeof(MohrCoulomb) );
+	self = (MohrCoulomb*) _YieldRheology_New( 
+			sizeOfSelf,
+			type,
+			_delete,
+			_print,
+			_copy,
+			_defaultConstructor,
+			_construct,
+			_build,
+			_initialise,
+			_execute,
+			_destroy,
+			_modifyConstitutiveMatrix,
+			_getYieldCriterion,
+			_getYieldIndicator,
+			_hasYielded,
+			name );
+	
+	/* Function pointers for this class that are not on the parent class should be set here */
+	
+	return self;
+}
+
+void _MohrCoulomb_Init(
+		MohrCoulomb*                        self,
+		FeVariable*                                        pressureField,
+		FeVariable*                                        velocityGradientsField,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		double                                             cohesion,
+		double                                             cohesionAfterSoftening,
+		double                                             frictionCoefficient,
+		double                                             frictionCoefficientAfterSoftening,
+		double                                             minimumYieldStress)
+{
+	StandardParticle                    materialPoint;
+	Dimension_Index                     dim   = materialPointsSwarm->dim;
+	
+	self->materialPointsSwarm     = materialPointsSwarm;
+	self->pressureField           = pressureField;
+	self->velocityGradientsField  = velocityGradientsField;
+	
+	self->cohesion = cohesion;
+	self->frictionCoefficient = frictionCoefficient;
+	
+	/* Strain softening of Cohesion - (linear weakening is assumed) */
+	/* needs a softening factor between +0 and 1 and a reference strain > 0 */
+	self->cohesionAfterSoftening = cohesionAfterSoftening;
+	
+	/* Strain softening of Friction - (linear weakening is assumed) */
+	/* needs a softening factor between +0 and 1 and a reference strain > 0 */
+	self->frictionCoefficientAfterSoftening = frictionCoefficientAfterSoftening;
+
+	self->minimumYieldStress = minimumYieldStress;
+}
+
+void* _MohrCoulomb_DefaultNew( Name name ) {
+	return (void*) _MohrCoulomb_New(
+		sizeof(MohrCoulomb),
+		MohrCoulomb_Type,
+		_YieldRheology_Delete,
+		_YieldRheology_Print,
+		_YieldRheology_Copy,
+		_MohrCoulomb_DefaultNew,
+		_MohrCoulomb_Construct,
+		_MohrCoulomb_Build,
+		_MohrCoulomb_Initialise,
+		_YieldRheology_Execute,
+		_YieldRheology_Destroy,
+		_MohrCoulomb_ModifyConstitutiveMatrix,
+		_MohrCoulomb_GetYieldCriterion,
+		_MohrCoulomb_GetYieldIndicator,
+		_MohrCoulomb_HasYielded,
+		name );
+}
+
+void _MohrCoulomb_Construct( void* rheology, Stg_ComponentFactory* cf ){
+	MohrCoulomb*   self           = (MohrCoulomb*)rheology;
+	FeVariable*                   pressureField;
+	MaterialPointsSwarm*          materialPointsSwarm;
+	FeVariable*                   velocityGradientsField;
+	
+	/* Construct Parent */
+	_YieldRheology_Construct( self, cf );
+
+	/* Make sure that there is strain weakening */
+	Journal_Firewall(
+		self->strainWeakening != NULL,
+		Journal_Register( Error_Type, self->type ),
+		"Error in func '%s' for %s '%s': MohrCoulomb rheology needs strain weakening.\n", 
+		__func__, self->type, self->name );
+	
+	materialPointsSwarm    = Stg_ComponentFactory_ConstructByKey( cf, self->name, 
+			"MaterialPointsSwarm", MaterialPointsSwarm, True );
+	pressureField          = Stg_ComponentFactory_ConstructByKey( cf, self->name,
+			"PressureField", FeVariable, True );
+	velocityGradientsField = Stg_ComponentFactory_ConstructByKey( cf, self->name,
+			"VelocityGradientsField", FeVariable, True );
+	
+	_MohrCoulomb_Init( 
+			self,
+			pressureField,
+			velocityGradientsField,
+			materialPointsSwarm, 
+			Stg_ComponentFactory_GetDouble( cf, self->name, "cohesion", 0.0 ),
+			Stg_ComponentFactory_GetDouble( cf, self->name, "cohesionAfterSoftening", 0.0 ),
+			Stg_ComponentFactory_GetDouble( cf, self->name, "frictionCoefficient", 0.0 ),
+			Stg_ComponentFactory_GetDouble( cf, self->name, "frictionCoefficientAfterSoftening", 0.0 ),
+			Stg_ComponentFactory_GetDouble( cf, self->name, "minimumYieldStress", 0.0 ) );
+}
+
+void _MohrCoulomb_Build( void* rheology, void* data ) {
+	MohrCoulomb*  self               = (MohrCoulomb*) rheology;
+
+	/* Build parent */
+	_YieldRheology_Build( self, data );
+}
+
+void _MohrCoulomb_Initialise( void* rheology, void* data ) {
+	MohrCoulomb*     self                  = (MohrCoulomb*) rheology;
+	double*                         ptr;
+	Particle_Index                  lParticle_I;
+	Particle_Index                  particleLocalCount;
+	double                          normalLength2;
+	double                          invNormalLength;
+	double                          initialDamageFraction;
+	MaterialPoint*                  materialPoint;
+	Index                           dof_I;
+	Dimension_Index                 dim                   = self->materialPointsSwarm->dim;
+	AbstractContext*                context = (AbstractContext*)data;
+
+	_YieldRheology_Initialise( self, data );
+	
+	/* Initialise these components now just in case this function is called before their own _Initialise function */
+	Stg_Component_Initialise( self->materialPointsSwarm, data, False );
+	Stg_Component_Initialise( self->strainWeakening, data, False );
+
+	/* We don't need to Initialise hasYieldedVariable because it's a parent variable and _YieldRheology_Initialise
+	 * has already been called */
+	particleLocalCount = self->hasYieldedVariable->variable->arraySize;
+
+	/* If restarting from checkpoint, don't change the parameters on the particles */
+	if ( !(context && (True == context->loadFromCheckPoint) ) ) {
+		for ( lParticle_I = 0 ; lParticle_I < particleLocalCount ; lParticle_I++ ) { 
+			Variable_SetValueChar( self->hasYieldedVariable->variable, lParticle_I, False );
+                }
+	}
+}
+	
+void _MohrCoulomb_ModifyConstitutiveMatrix( 
+		void*                                              rheology, 
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi )
+{
+	MohrCoulomb*     self                  = (MohrCoulomb*) rheology;
+	double                          postFailureWeakening;
+	double                          yieldCriterion;
+	double                          yieldIndicator;  /* A materialPoint will yield if yieldCriterion < yieldIndicator */
+
+	/* Don't want to yield on the first ever solve */
+	if ( constitutiveMatrix->previousSolutionExists == False ) {
+		return;
+	}
+
+	/* Calculate and store values of stress, pressure, velocity gradients, eigenvectors so they are only calculated once */
+	_MohrCoulomb_StoreCurrentParameters( self, constitutiveMatrix, materialPointsSwarm, lElement_I, materialPoint, xi );
+
+	_YieldRheology_ModifyConstitutiveMatrix( self, constitutiveMatrix, materialPointsSwarm, lElement_I, materialPoint, xi );
+}
+
+double _MohrCoulomb_GetYieldIndicator( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi )
+{
+	MohrCoulomb*           self             = (MohrCoulomb*) rheology;
+	double*                               stress           = self->currentStress;
+	Eigenvector*                          eigenvectorList  = self->currentEigenvectorList;
+	double                                sigma_ns;
+	Dimension_Index                       dim              = constitutiveMatrix->dim;
+	double                                stressMin;
+	double                                stressMax;
+	double                                theta;
+
+        /* Failure criterion from stress field - we can calculate
+         * whether the material has failed at the current stress from
+         * the principle stresses */
+        theta = 0.5 * atan( 1.0/ _MohrCoulomb_EffectiveFrictionCoefficient( self, materialPoint ) );
+        
+        stressMin = eigenvectorList[0].eigenvalue;
+        stressMax = (dim == 2 ? eigenvectorList[1].eigenvalue : eigenvectorList[2].eigenvalue);
+        
+        sigma_ns = 0.5 *  sin( 2.0 * theta ) * ( stressMax - stressMin );
+	
+	return fabs(sigma_ns);
+}
+
+double _MohrCoulomb_GetYieldCriterion( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi )
+{
+	MohrCoulomb*          self             = (MohrCoulomb*) rheology;
+	double                               minimumYieldStress;
+	double                               effectiveCohesion;
+	double                               effectiveFrictionCoefficient;
+	double                               frictionalStrength;
+	double                               sigma_nn;
+		
+	/* Calculate frictional strength */
+	effectiveFrictionCoefficient = _MohrCoulomb_EffectiveFrictionCoefficient( self, materialPoint );
+	effectiveCohesion            = _MohrCoulomb_EffectiveCohesion( self, materialPoint );
+	sigma_nn                     = _MohrCoulomb_Sigma_nn( self, materialPoint, constitutiveMatrix->dim );
+
+	frictionalStrength = effectiveFrictionCoefficient * sigma_nn + effectiveCohesion;
+
+	/* Make sure frictionalStrength is above the minimum */
+	minimumYieldStress = self->minimumYieldStress;
+	
+	if ( frictionalStrength < minimumYieldStress) 
+		frictionalStrength = minimumYieldStress;
+	
+	return frictionalStrength;
+}
+
+/* This sets the viscosity to yield_stress/(2*second stress invariant) if
+   the material has failed. */
+
+void _MohrCoulomb_HasYielded( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		double                                             yieldCriterion,
+		double                                             yieldIndicator )
+{
+  MohrCoulomb*                     self = (MohrCoulomb*)rheology;
+  Eigenvector*                     eigenvectorList  = self->currentEigenvectorList;
+  Dimension_Index                  dim = constitutiveMatrix->dim;
+  double viscosity;
+  
+  viscosity = yieldCriterion/(2*self->strainRateSecondInvariant);
+
+  ConstitutiveMatrix_SetIsotropicViscosity( constitutiveMatrix, viscosity );
+}
+
+double _MohrCoulomb_EffectiveCohesion( void* rheology, void* materialPoint ) {
+	MohrCoulomb*    self                          = (MohrCoulomb*) rheology;
+	double                         effectiveCohesion;
+	
+        double strainWeakeningRatio = StrainWeakening_CalcRatio( self->strainWeakening, materialPoint );
+	
+        effectiveCohesion =  self->cohesion * (1.0 - strainWeakeningRatio) + 
+          self->cohesionAfterSoftening * strainWeakeningRatio;
+
+	return effectiveCohesion;
+}
+
+double _MohrCoulomb_EffectiveFrictionCoefficient( void* rheology, void* materialPoint ) {
+	MohrCoulomb*    self                          = (MohrCoulomb*) rheology;
+	static double                  effectiveFrictionCoefficient  = 0.0;
+
+        double strainWeakeningRatio = StrainWeakening_CalcRatio( self->strainWeakening, materialPoint );
+	
+        effectiveFrictionCoefficient = self->frictionCoefficient * (1.0 - strainWeakeningRatio) +
+          self->frictionCoefficientAfterSoftening * strainWeakeningRatio;	
+
+	return effectiveFrictionCoefficient;
+}
+
+double _MohrCoulomb_Sigma_nn( void* rheology, void* materialPoint, Dimension_Index dim ) {
+	MohrCoulomb*      self             = (MohrCoulomb*) rheology;
+	double                           pressure         = self->currentPressure;
+	Eigenvector*                     eigenvectorList  = self->currentEigenvectorList;
+	double                           tau_nn;
+	double                           sigma_nn;
+
+        double stressMin = eigenvectorList[0].eigenvalue;
+        double stressMax = (dim == 2 ? eigenvectorList[1].eigenvalue : eigenvectorList[2].eigenvalue);
+        double theta = 0.5 * atan( 1.0/ _MohrCoulomb_EffectiveFrictionCoefficient( self, materialPoint ) );
+        
+        tau_nn = 0.5 * (( stressMax + stressMin ) + cos( 2.0 * theta ) * ( stressMax - stressMin ));
+
+	sigma_nn = pressure - tau_nn;
+
+	return sigma_nn;
+}
+
+void _MohrCoulomb_StoreCurrentParameters( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix, 
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I, 
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi ) 
+{
+	MohrCoulomb*          self               = (MohrCoulomb*) rheology;
+	SymmetricTensor                      strainRate;
+	Dimension_Index                      dim                = constitutiveMatrix->dim;
+        Eigenvector                          evectors[3];
+	
+	FeVariable_InterpolateWithinElement( self->pressureField, lElement_I, xi, &self->currentPressure );	
+	FeVariable_InterpolateWithinElement( self->velocityGradientsField, lElement_I, xi, self->currentVelocityGradient );
+	
+	TensorArray_GetSymmetricPart( self->currentVelocityGradient, dim, strainRate );
+	ConstitutiveMatrix_CalculateStress( constitutiveMatrix, strainRate, self->currentStress );
+	
+	SymmetricTensor_CalcAllEigenvectors( self->currentStress, dim, self->currentEigenvectorList );
+
+        SymmetricTensor_CalcAllEigenvectors( strainRate, dim, evectors);
+        if(dim == 2)
+          {
+            self->strainRateSecondInvariant=0.5*fabs(evectors[0].eigenvalue-evectors[1].eigenvalue);
+          }
+}


Property changes on: long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.c
___________________________________________________________________
Name: svn:keywords
   + LastChangedDate Author Id
Name: svn:eol-style
   + native

Added: long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.h
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.h	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.h	2006-10-11 18:39:23 UTC (rev 4787)
@@ -0,0 +1,152 @@
+/*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+** Copyright (c) 2005, Monash Cluster Computing 
+** All rights reserved.
+** Redistribution and use in source and binary forms, with or without modification,
+** are permitted provided that the following conditions are met:
+**
+** 		* Redistributions of source code must retain the above copyright notice, 
+** 			this list of conditions and the following disclaimer.
+** 		* Redistributions in binary form must reproduce the above copyright 
+**			notice, this list of conditions and the following disclaimer in the 
+**			documentation and/or other materials provided with the distribution.
+** 		* Neither the name of the Monash University nor the names of its contributors 
+**			may be used to endorse or promote products derived from this software 
+**			without specific prior written permission.
+**
+** THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 
+** "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, 
+** THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR 
+** PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS 
+** BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR 
+** CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF 
+** SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) 
+** HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 
+** LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT 
+** OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
+**
+**
+** Contact:
+*%		Louis Moresi - Louis.Moresi at sci.monash.edu.au
+*%
+** Contributors:
+*+		Robert Turnbull
+*+		Vincent Lemiale
+*+		Louis Moresi
+*+		David May
+*+		David Stegman
+*+		Mirko Velic
+*+		Patrick Sunter
+*+		Julian Giordani
+*+
+** $Id$
+** 
+**~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
+
+
+#ifndef __Underworld_Rheology_MohrCoulomb_h__
+#define __Underworld_Rheology_MohrCoulomb_h__
+
+	/** Textual name of this class - This is a global pointer which is used for times when you need to refer to class and not a particular instance of a class */
+	extern const Type MohrCoulomb_Type;
+	
+	/** Rheology class contents - this is defined as a macro so that sub-classes of this class can use this macro at the start of the definition of their struct */
+	#define __MohrCoulomb \
+		/* Parent info */ \
+ 		__YieldRheology \
+		/* Virtual functions go here */ \
+		/* Material Parameters */\
+		ExtensionInfo_Index                                 particleExtHandle;                     \
+		/* Param passed in */ \
+		MaterialPointsSwarm*                                materialPointsSwarm;                                 \
+		FeVariable*                                         pressureField;                         \
+		FeVariable*                                         velocityGradientsField;                \
+		/* Director component is used to update the normal */\
+		Director*                                           director;                              \
+		double                                              cohesion;                              \
+		double                                              cohesionAfterSoftening;                \
+		double                                              frictionCoefficient;                   \
+		double                                              frictionCoefficientAfterSoftening;     \
+		double                                              minimumYieldStress;                    \
+		/* Stored values that are calculated once for each particle */ \
+		Eigenvector                                         currentEigenvectorList[3];             \
+		TensorArray                                         currentVelocityGradient;               \
+		SymmetricTensor                                     currentStress;                         \
+		double                                              currentPressure;                       \
+		double                                              strainRateSecondInvariant;                   \
+
+	
+	struct MohrCoulomb { __MohrCoulomb };
+	
+	/** Private Constructor: This will accept all the virtual functions for this class as arguments. */
+	MohrCoulomb* _MohrCoulomb_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,
+		Rheology_ModifyConstitutiveMatrixFunction*         _modifyConstitutiveMatrix,
+		YieldRheology_GetYieldCriterionFunction*           _getYieldCriterion,
+		YieldRheology_GetYieldIndicatorFunction*           _getYieldIndicator,
+		YieldRheology_HasYieldedFunction*                  _hasYielded,
+		Name                                               name ) ;
+	
+	/* 'Stg_Component' implementations */
+	void* _MohrCoulomb_DefaultNew( Name name ) ;
+	void _MohrCoulomb_Construct( void* rheology, Stg_ComponentFactory* cf );
+
+	void _MohrCoulomb_Build( void* rheology, void* data );
+	void _MohrCoulomb_Initialise( void* rheology, void* data ) ;
+
+	/* 'YieldRheology' implementations */
+	void _MohrCoulomb_ModifyConstitutiveMatrix( 
+		void*                                              rheology, 
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi );
+		
+	double _MohrCoulomb_GetYieldCriterion( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi );
+	
+	double _MohrCoulomb_GetYieldIndicator( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi );
+	
+	void _MohrCoulomb_HasYielded( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix,
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I,
+		MaterialPoint*                                     materialPoint,
+		double                                             yieldCriterion,
+		double                                             yieldIndicator );
+	
+	double _MohrCoulomb_EffectiveCohesion( void* rheology, void* materialPoint ) ;
+	double _MohrCoulomb_EffectiveFrictionCoefficient( void* rheology, void* materialPoint ) ;
+	double _MohrCoulomb_Sigma_nn( void* rheology, void* materialPoint, Dimension_Index dim ) ;
+
+	void _MohrCoulomb_StoreCurrentParameters( 
+		void*                                              rheology,
+		ConstitutiveMatrix*                                constitutiveMatrix, 
+		MaterialPointsSwarm*                               materialPointsSwarm,
+		Element_LocalIndex                                 lElement_I, 
+		MaterialPoint*                                     materialPoint,
+		Coord                                              xi ) ;
+	
+#endif


Property changes on: long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.h
___________________________________________________________________
Name: svn:keywords
   + LastChangedDate Author Id
Name: svn:eol-style
   + native

Added: long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.meta
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.meta	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/MohrCoulomb.meta	2006-10-11 18:39:23 UTC (rev 4787)
@@ -0,0 +1,84 @@
+<?xml version="1.0"?>
+<!DOCTYPE StGermainData SYSTEM "stgermain.dtd">
+<StGermainData xmlns="http://www.vpac.org/StGermain/XML_IO_Handler/Jun2003">
+
+<param name="Name">MohrCoulomb</param>
+<param name="Organisation">MCC</param>
+<param name="Project">Underworld</param>
+<param name="Location">./Underworld/Rheology/src/</param>
+<param name="Project Web">http://mcc.monash.edu.au/Underworld</param>
+<param name="Copyright">Copyright (c) 2005, Monash Cluster Computing</param>
+<param name="License">http://www.opensource.org/licenses/bsd-license.php</param>
+<param name="Parent">YieldRheology</param>
+<param name="Description">Implements faulting behaviour based on Mohr-Coulomb failure criterion. See following reference for details : Anisotropic viscous models of large-deformation Mohr–Coulomb failure, Moresi, L.; Mühlhaus, H.-B., Philosophical Magazine, Volume 86, Numbers 21-22, -22/21 July–1 August 2006, pp. 3287-3305</param>
+
+<!--Now the interesting stuff-->
+
+<list name="Params">
+	<struct>
+		<param name="Name">cohesion</param>
+		<param name="Type">Double</param>
+		<param name="Default">0.0</param>
+		<param name="Description">Cohesion for a pristine material</param>
+	</struct>
+	<struct>
+		<param name="Name">cohesionAfterSoftening</param>
+		<param name="Type">Double</param>
+		<param name="Default">0.0</param>
+		<param name="Description">Value of the cohesion after the material has weakened</param>
+	</struct>
+	<struct>
+		<param name="Name">frictionCoefficient</param>
+		<param name="Type">Double</param>
+		<param name="Default">0.0</param>
+		<param name="Description">Defines the dependance of the criterion to the normal stress</param>
+	</struct>
+	<struct>
+		<param name="Name">frictionCoefficientAfterSoftening</param>
+		<param name="Type">Double</param>
+		<param name="Default">0.0</param>
+		<param name="Description">This is the value of the frictionCoefficient after the material has weakened</param>
+	</struct>
+	<struct>
+		<param name="Name">minimumYieldStress</param>
+		<param name="Type">Double</param>
+		<param name="Default">0.0</param>
+		<param name="Description">The value of the failure criterion (defining the yield surface) can not be smaller than this parameter</param>
+	</struct>
+</list>
+
+<list name="Dependencies">
+
+	<struct>
+		<param name="Essential">True</param>
+		<param name="Type">FeVariable</param>
+		<param name="Fallback Key">PressureField</param>
+		<param name="Description">Essential since the MohrCoulomb failure criterion is pressure dependant</param>
+	</struct>
+
+	<struct>
+		<param name="Essential">True</param>
+		<param name="Type">IntegrationPointsSwarm</param>
+		<param name="Description">Defines the set of material points</param>
+	</struct>
+	
+	<struct>
+		<param name="Essential">True</param>
+		<param name="Type">FeVariable</param>
+		<param name="Fallback Key">VelocityGradientsFieldField</param>
+		<param name="Description">The resolved component of this tensor wrt the failure plane are mainly used to determine the most favourable oriented plane for failure to occur </param>
+	</struct>
+
+	<struct>
+		<param name="Essential">True</param>
+		<param name="Type">StrainWeakening</param>
+		<param name="Description">Define the weakening behaviour of material parameters. This component dependency is first defined in YieldRheology class, but as non essential. Since we want it for MC criterion, we make it essential here.</param>
+	</struct>
+
+
+</list>
+<!-- Add an example XML if possible -->
+<param name="Example">Underworld/InputFiles/ExtensionFMM.xml</param>
+<param name="Example">Underworld/InputFiles/ExtensionFMM3D.xml</param>
+
+</StGermainData>

Modified: long/3D/Gale/trunk/src/Underworld/Rheology/src/Rheology.h
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/Rheology.h	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/Rheology.h	2006-10-11 18:39:23 UTC (rev 4787)
@@ -67,6 +67,7 @@
 	#include "Byerlee.h"
 	#include "DruckerPrager.h"
 	#include "FaultingMoresiMuhlhaus2006.h"
+	#include "MohrCoulomb.h"
 	#include "StrainWeakening.h"
 	#include "BuiterStrainWeakening.h"
 	#include "Director.h"

Modified: long/3D/Gale/trunk/src/Underworld/Rheology/src/SConscript
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/SConscript	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/SConscript	2006-10-11 18:39:23 UTC (rev 4787)
@@ -35,6 +35,7 @@
 FrankKamenetskii.h
 Init.h
 MaterialViscosity.h
+MohrCoulomb.h
 MultiRheologyMaterial.h
 NonNewtonian.h
 OrthotropicAligned.h
@@ -67,6 +68,7 @@
 FrankKamenetskii.c
 Init.c
 MaterialViscosity.c
+MohrCoulomb.c
 MultiRheologyMaterial.c
 NonNewtonian.c
 OrthotropicAligned.c
@@ -95,6 +97,7 @@
 FaultingMoresiMuhlhaus2006.meta
 FrankKamenetskii.meta
 MaterialViscosity.meta
+MohrCoulomb.meta
 MultiRheologyMaterial.meta
 NonNewtonian.meta
 OrthotropicAligned.meta

Modified: long/3D/Gale/trunk/src/Underworld/Rheology/src/types.h
===================================================================
--- long/3D/Gale/trunk/src/Underworld/Rheology/src/types.h	2006-10-11 18:39:16 UTC (rev 4786)
+++ long/3D/Gale/trunk/src/Underworld/Rheology/src/types.h	2006-10-11 18:39:23 UTC (rev 4787)
@@ -68,8 +68,9 @@
 	typedef struct Byerlee                          Byerlee;
 	
 	typedef struct DruckerPrager                    DruckerPrager;
-	typedef struct FaultingMoresiMuhlhaus2006       FaultingMoresiMuhlhaus2006;
-	
+	typedef struct FaultingMoresiMuhlhaus2006       FaultingMoresiMuhlhaus2006;	
+	typedef struct MohrCoulomb                      MohrCoulomb;
+
 	typedef struct StrainWeakening                  StrainWeakening;
 	typedef struct BuiterStrainWeakening		BuiterStrainWeakening;
 	typedef struct Director                         Director;



More information about the cig-commits mailing list