SPModel: plugins/isostaticFlexure/IsostaticFlexure.c Source File
VPAC - Computational Software Development
Main | SPModel | StGermain FrameWork |
Main Page | Alphabetical List | Class List | Directories | File List | Class Members | File Members

IsostaticFlexure.c

Go to the documentation of this file.
00001 /*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00002 **
00003 ** Copyright (C), 2004, Victorian Partnership for Advanced Computing (VPAC) Ltd, 110 Victoria Street, Melbourne, 3053, Australia.
00004 **
00005 ** Authors:
00006 **  Ogar R. Widjaja, Computational Scientist, VPAC.
00007 **  Raquibul Hassan, Software Engineer, VPAC. (raq@vpac.org)
00008 **  Keith Hsuan, Computational Scientist, VPAC (keith@vpac.org)
00009 **  William F. Appelbe, Director, VPAC. (bill@vpac.org)
00010 **  Stevan M. Quenette, Senior Software Engineer, VPAC. (steve@vpac.org)
00011 **  Patrick D. Sunter, Software Engineer, VPAC. (patrick@vpac.org)
00012 **
00013 ** This file may be distributed under the terms of the VPAC Public License
00014 ** as defined by VPAC of Australia and appearing in the file
00015 ** LICENSE.VPL included in the packaging of this file.
00016 **
00017 ** This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
00018 ** WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
00019 **
00020 ** $Id: IsostaticFlexure.c 225 2005-12-22 00:01:19Z AlanLo $
00021 **
00022 **~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~*/
00023 
00024 
00025 #include <mpi.h>
00026 #include <StGermain/StGermain.h>
00027 #include "SPModel/SPModel.h"
00028 #include "flex2d.h"
00029 #include <stdlib.h>
00030 #include <string.h>
00031 #include <assert.h>
00032 #include "IsostaticFlexure.h"
00033 
00034 /* Textual name of this module */
00035 const Type SPModelIsostaticFlexure_Type = "SPModelIsostaticFlexure";
00036 
00037 ExtensionInfo_Index SPModel_Simulation_ContextExtHandle;
00038 
00039 void* _SPModelIsostaticFlexure_DefaultNew( Name name ) {
00040     return Codelet_New(
00041             SPModelIsostaticFlexure_Type,
00042             _SPModelIsostaticFlexure_DefaultNew,
00043             _SPModelIsostaticFlexure_Construct,
00044             _Codelet_Build,
00045             _Codelet_Initialise,
00046             _Codelet_Execute,
00047             _Codelet_Destroy,
00048             name );
00049 }
00050 
00051 void _SPModelIsostaticFlexure_Construct( void* component, Stg_ComponentFactory* data ) {    
00052     SPModel_Context*                context;
00053 
00054     context = (SPModel_Context*)Stg_ComponentFactory_ConstructByName( data, "context", SPModel_Context, True );
00055     
00056     assert( context );
00057 
00058     SPModel_Simulation_ContextExtHandle = ExtensionManager_GetHandle( context->extensionMgr, "Simulation" );
00059 
00060     _SPModelIsostaticFlexure_Init( context );
00061     
00062     Journal_DPrintf( context->debug, "In %s():\n", __func__ );
00063     ContextEP_Append( context, AbstractContext_EP_Solve, SPModelIsostaticFlexure_Solve );
00064 }   
00065 
00066 void _SPModelIsostaticFlexure_Init( void* _context ) {
00067 
00068 }
00069 
00070 float interpolateFunction ( int meshLocalIndex, void *args )
00071 {
00072     float result = 0.0f;
00073     SPModel_Context *context = (SPModel_Context *) args;
00074     SPModelSimulationContextExtension *simulationExt = NULL;
00075 
00076     assert( context );
00077 
00078     simulationExt = (SPModelSimulationContextExtension*)ExtensionManager_Get
00079                 ( context->extensionMgr, context, SPModel_Simulation_ContextExtHandle );
00080 
00081     assert( simulationExt );
00082 
00083     result = 0.0f;
00084     
00085         result = -simulationExt->erosion[meshLocalIndex] * context->bedRockDensity*context->gravity
00086                  +simulationExt->sedimentHistory[meshLocalIndex] * context->sedimentDensity*context->gravity
00087                  +simulationExt->uplift[meshLocalIndex] * context->bedRockDensity*context->gravity;
00088     
00089     return result;
00090 }
00091 
00092 void SPModelIsostaticFlexure_Solve( void *_context ) {
00093     double **loadg;
00094     double **deflec;
00095     double ***loadt;
00096     double **speq;
00097     SPModelSimulationContextExtension *simulationExt = NULL;
00098     int i, j;
00099 
00100     SPModel_Context *context = NULL;
00101     SurfaceMesh *globalMesh = NULL, *localMesh = NULL;
00102     SurfaceRegularMesh *regularMesh = NULL;
00103     _Interpolator *interpolator = NULL;
00104     
00105     context = (SPModel_Context*)_context;
00106     assert( context );
00107     
00108     if( ( context->timeStep % context->flexureInterval ) ){
00109         return;
00110     }
00111 
00112     simulationExt = (SPModelSimulationContextExtension*)ExtensionManager_Get
00113                     ( context->extensionMgr, context, SPModel_Simulation_ContextExtHandle );
00114     
00115     globalMesh = context->globalMesh;
00116     localMesh = context->localMesh;
00117     interpolator = context->interpolator;
00118     regularMesh = context->regularMesh;
00119     
00120     _Interpolator_InterpolateFromMeshToGrid( interpolator, interpolateFunction, context );
00121     
00122     /*Stg_Class_Print( interpolator, Journal_Register( Info_Type, "lota" ) );*/
00123 
00124     SurfaceRegularMesh_GatherData( regularMesh, interpolator );
00125 
00126     if( localMesh->rank == MASTER_PROC ){
00127 
00128         loadg = dmatrix(1,regularMesh->ny,1,regularMesh->nx);
00129         deflec = dmatrix(1,regularMesh->ny,1,regularMesh->nx);
00130         loadt = d3tensor(1,1,1,2*regularMesh->nx,1,2*regularMesh->ny);
00131         speq = dmatrix(1,1,1,2*2*regularMesh->nx);
00132 
00133         for( i=1; i<=regularMesh->nx; i++ ){
00134             for( j=1; j<=regularMesh->ny; j++ ){
00135                 loadg[i][j] = regularMesh->dataArray[i-1][j-1];
00136             }
00137         }
00138 
00139         cosfilt( context->elasticPlateLength, regularMesh->nx, context->elasticPlateLength, 
00140                 regularMesh->ny, loadg, context->flexuralRigidity,
00141                 context->asthenosphereDensity, context->gravity, deflec, loadt, speq );
00142 
00143         for( i=1; i<=regularMesh->nx; i++ ){
00144             for( j=1; j<=regularMesh->nx; j++ ){
00145                 regularMesh->dataArray[i-1][j-1] = deflec[i][j];
00146             }
00147         }
00148     
00149         free_dmatrix(loadg,1,regularMesh->ny,1,regularMesh->nx);
00150         free_dmatrix(deflec,1,regularMesh->ny,1,regularMesh->nx);
00151         free_d3tensor(loadt,1,1,1,2*regularMesh->nx,1,2*regularMesh->ny);
00152         free_dmatrix(speq,1,1,1,2*2*regularMesh->nx);
00153     }
00154     
00155     SurfaceRegularMesh_BroadcastData( regularMesh );
00156     
00157     _Interpolator_InterpolateFromGridToMesh( interpolator, regularMesh->dataArray, simulationExt->hiso );
00158 
00159     for( i=0; i<localMesh->numNodes; i++ ){
00160         simulationExt->hisotot[i] = -simulationExt->hisoPrev[i];
00161         simulationExt->hisotot[i] += simulationExt->hiso[i];
00162         simulationExt->hisoPrev[i] = simulationExt->hiso[i];
00163     }
00164 
00165     for( i=0; i<localMesh->myLoad; i++ ){
00166         localMesh->h[i] += simulationExt->hisotot[localMesh->id[i]-1];
00167     }
00168 }
00169 
00170 Index SPModelIsostaticFlexure_Register( PluginsManager* pluginsManager ) {
00171     return PluginsManager_Submit( pluginsManager, SPModelIsostaticFlexure_Type, "0", _SPModelIsostaticFlexure_DefaultNew );
00172 }
00173