/MultiPhysics/Poroelasticity/FiniteElasticityDarcy/IncompressibleElasticityDrivenDarcyIO/src/IncompressibleElasticityDrivenDarcyIOExample.f90
FORTRAN Modern | 1137 lines | 571 code | 184 blank | 382 comment | 0 complexity | fa0189050314fecf43370947ff8449ab MD5 | raw file
Large files files are truncated, but you can click here to view the full file
1!> \file 2!> \author Christian Michler, Adam Reeve 3!> \brief This is an example program to solve a coupled Finite Elastiticity Darcy equation using openCMISS calls. 4!> 5!> \section LICENSE 6!> 7!> Version: MPL 1.1/GPL 2.0/LGPL 2.1 8!> 9!> The contents of this file are subject to the Mozilla Public License 10!> Version 1.1 (the "License"); you may not use this file except in 11!> compliance with the License. You may obtain a copy of the License at 12!> http://www.mozilla.org/MPL/ 13!> 14!> Software distributed under the License is distributed on an "AS IS" 15!> basis, WITHOUT WARRANTY OF ANY KIND, either express or implied. See the 16!> License for the specific language governing rights and limitations 17!> under the License. 18!> 19!> The Original Code is OpenCMISS 20!> 21!> The Initial Developer of the Original Code is University of Auckland, 22!> Auckland, New Zealand and University of Oxford, Oxford, United 23!> Kingdom. Portions created by the University of Auckland and University 24!> of Oxford are Copyright (C) 2007 by the University of Auckland and 25!> the University of Oxford. All Rights Reserved. 26!> 27!> Contributor(s): 28!> 29!> Alternatively, the contents of this file may be used under the terms of 30!> either the GNU General Public License Version 2 or later (the "GPL"), or 31!> the GNU Lesser General Public License Version 2.1 or later (the "LGPL"), 32!> in which case the provisions of the GPL or the LGPL are applicable instead 33!> of those above. If you wish to allow use of your version of this file only 34!> under the terms of either the GPL or the LGPL, and not to allow others to 35!> use your version of this file under the terms of the MPL, indicate your 36!> decision by deleting the provisions above and replace them with the notice 37!> and other provisions required by the GPL or the LGPL. If you do not delete 38!> the provisions above, a recipient may use your version of this file under 39!> the terms of any one of the MPL, the GPL or the LGPL. 40!> 41 42!> \example MultiPhysics/Poroelasticity/FiniteElasticityDarcy/IncompressibleElasticityDrivenDarcy/src/IncompressibleElasticityDrivenDarcyExample.f90 43!! Example program to solve coupled FiniteElasticityDarcy equations using OpenCMISS calls. 44!! \par Latest Builds: 45!! \li <a href='http://autotest.bioeng.auckland.ac.nz/opencmiss-build/logs_x86_64-linux/MultiPhysics/Poroelasticity/FiniteElasticityDarcy/IncompressibleElasticityDrivenDarcy/build-intel'>Linux Intel Build</a> 46!! \li <a href='http://autotest.bioeng.auckland.ac.nz/opencmiss-build/logs_x86_64-linux/MultiPhysics/Poroelasticity/FiniteElasticityDarcy/IncompressibleElasticityDrivenDarcy/build-intel'>Linux GNU Build</a> 47!! 48!< 49 50! ! 51! ! This example considers a coupled Finite Elasticity Darcy problem 52! ! 53 54!> Main program 55 56PROGRAM FINITEELASTICITYDARCYIOEXAMPLE 57 58 ! 59 !================================================================================================================================ 60 ! 61 62 !PROGRAM LIBRARIES 63 64 USE OPENCMISS 65! USE FLUID_MECHANICS_IO_ROUTINES 66 USE IOSTUFF 67 USE MPI 68 69#ifdef WIN32 70 USE IFQWINCMISS 71#endif 72 73 ! 74 !================================================================================================================================ 75 ! 76 77 !PROGRAM VARIABLES AND TYPES 78 79 IMPLICIT NONE 80 81 !Test program parameters 82 83 REAL(CMISSDP), PARAMETER :: Y_DIM=1.0_CMISSDP 84 REAL(CMISSDP), PARAMETER :: X_DIM=1.0_CMISSDP 85 REAL(CMISSDP), PARAMETER :: Z_DIM=1.0_CMISSDP 86 87 INTEGER(CMISSIntg), PARAMETER :: LinearBasisUserNumber=1 88 INTEGER(CMISSIntg), PARAMETER :: QuadraticBasisUserNumber=2 89 INTEGER(CMISSIntg), PARAMETER :: CubicBasisUserNumber=3 90 91 INTEGER(CMISSIntg), PARAMETER :: CoordinateSystemUserNumber=1 92 INTEGER(CMISSIntg), PARAMETER :: RegionUserNumber=2 93 INTEGER(CMISSIntg), PARAMETER :: MeshUserNumber=3 94 INTEGER(CMISSIntg), PARAMETER :: DecompositionUserNumber=4 95 INTEGER(CMISSIntg), PARAMETER :: GeometricFieldUserNumber=5 96 INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDarcy=8 97 INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumberDarcy=12 98 INTEGER(CMISSIntg), PARAMETER :: ProblemUserNumber=14 99 INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldUserNumberDarcy=22 100 101 INTEGER(CMISSIntg), PARAMETER :: ControlLoopSolidNumber=1 102 INTEGER(CMISSIntg), PARAMETER :: ControlLoopFluidNumber=2 103 INTEGER(CMISSIntg), PARAMETER :: ControlLoopSubiterationNumber=1 104 INTEGER(CMISSIntg), PARAMETER :: SolverSolidIndex=1 105 INTEGER(CMISSIntg), PARAMETER :: SolverDarcyIndex=1 106 INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDarcyPorosity=1 107 INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDarcyPermOverVis=2 108 109 INTEGER(CMISSIntg), PARAMETER :: FieldGeometryNumberOfVariables=1 110 INTEGER(CMISSIntg), PARAMETER :: FieldGeometryNumberOfComponents=3 111 112 !Program types 113 114 !Program variables 115 116 INTEGER(CMISSIntg) :: NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS,NUMBER_GLOBAL_Z_ELEMENTS 117 118! INTEGER(CMISSIntg) :: MPI_IERROR 119 INTEGER(CMISSIntg) :: NumberOfComputationalNodes,NumberOfDomains,ComputationalNodeNumber 120 121 INTEGER(CMISSIntg) :: NUMBER_OF_DIMENSIONS 122 123 INTEGER(CMISSIntg) :: MAXIMUM_ITERATIONS 124 INTEGER(CMISSIntg) :: RESTART_VALUE 125 126 INTEGER(CMISSIntg) :: EQUATIONS_DARCY_OUTPUT 127 INTEGER(CMISSIntg) :: COMPONENT_NUMBER 128 129 INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_DARCY_OUTPUT_FREQUENCY 130 INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_DARCY_OUTPUT_TYPE 131 INTEGER(CMISSIntg) :: LINEAR_SOLVER_DARCY_OUTPUT_TYPE 132 INTEGER(CMISSIntg) :: LINEAR_SOLVER_MAT_PROPERTIES_OUTPUT_TYPE 133 134 REAL(CMISSDP) :: GEOMETRY_TOLERANCE 135 INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_SOLID 136 REAL(CMISSDP) :: INITIAL_FIELD_DARCY(4) 137 REAL(CMISSDP) :: DIVERGENCE_TOLERANCE 138 REAL(CMISSDP) :: RELATIVE_TOLERANCE 139 REAL(CMISSDP) :: ABSOLUTE_TOLERANCE 140 REAL(CMISSDP) :: LINESEARCH_ALPHA 141 REAL(CMISSDP) :: VALUE 142 REAL(CMISSDP) :: POROSITY_PARAM_DARCY, PERM_OVER_VIS_PARAM_DARCY 143 144 LOGICAL :: EXPORT_FIELD_IO 145 LOGICAL :: LINEAR_SOLVER_DARCY_DIRECT_FLAG 146 147 !CMISS variables 148 149 !Regions 150 TYPE(CMISSRegionType) :: Region 151 TYPE(CMISSRegionType) :: WorldRegion 152 !Coordinate systems 153 TYPE(CMISSCoordinateSystemType) :: CoordinateSystem 154 TYPE(CMISSCoordinateSystemType) :: WorldCoordinateSystem 155 !Basis 156! TYPE(CMISSBasisType) :: CubicBasis, QuadraticBasis, LinearBasis, Bases(2) 157 TYPE(CMISSBasisType),allocatable :: Bases(:) 158 !Meshes 159 TYPE(CMISSMeshType) :: Mesh 160 TYPE(CMISSGeneratedMeshType) :: GeneratedMesh 161 162 !Decompositions 163 TYPE(CMISSDecompositionType) :: Decomposition 164 !Fields 165 TYPE(CMISSFieldsType) :: Fields 166 !Field types 167 TYPE(CMISSFieldType) :: GeometricField 168 TYPE(CMISSFieldType) :: MaterialsFieldDarcy 169 TYPE(CMISSFieldType) :: EquationsSetFieldDarcy 170 !Boundary conditions 171 TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsDarcy 172 !Equations sets 173 TYPE(CMISSEquationsSetType) :: EquationsSetDarcy 174 !Equations 175 TYPE(CMISSEquationsType) :: EquationsDarcy 176 !Problems 177 TYPE(CMISSProblemType) :: Problem 178 !Control loops 179 TYPE(CMISSControlLoopType) :: ControlLoop 180 !Solvers 181 TYPE(CMISSSolverType) :: DynamicSolverDarcy 182 TYPE(CMISSSolverType) :: LinearSolverDarcy 183! TYPE(CMISSSolverType) :: LinearSolverSolid 184 !Solver equations 185 TYPE(CMISSSolverEquationsType) :: SolverEquationsDarcy 186 ! nodes and elements 187 TYPE(CMISSNodesType) :: Nodes 188 TYPE(CMISSMeshElementsType),allocatable :: Elements(:) 189 190 !Other variables 191 INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face1Nodes(:),Face2Nodes(:) 192 INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face3Nodes(:),Face4Nodes(:) 193 INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face5Nodes(:),Face6Nodes(:) 194 INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face7Nodes(:),Face8Nodes(:) 195 INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face9Nodes(:),Face10Nodes(:) 196 INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face11Nodes(:),Face12Nodes(:) 197 INTEGER(CMISSIntg) :: FaceXi(6) 198 INTEGER(CMISSIntg) :: NN,NODE,NodeDomain 199 REAL(CMISSDP) :: XCoord,YCoord,ZCoord 200 LOGICAL :: X_FIXED,Y_FIXED !,X_OKAY,Y_OKAY 201 202#ifdef WIN32 203 !Quickwin type 204 LOGICAL :: QUICKWIN_STATUS=.FALSE. 205 TYPE(WINDOWCONFIG) :: QUICKWIN_WINDOW_CONFIG 206#endif 207 208 !Generic CMISS variables 209 210 INTEGER(CMISSIntg) :: EquationsSetIndex 211 INTEGER(CMISSIntg) :: Err 212 213 214 INTEGER(CMISSIntg) :: DIAG_LEVEL_LIST(5) 215! CHARACTER(LEN=255) :: DIAG_ROUTINE_LIST(8) !,TIMING_ROUTINE_LIST(1) 216 CHARACTER(LEN=255) :: DIAG_ROUTINE_LIST(1) !,TIMING_ROUTINE_LIST(1) 217 218 ! 219 !-------------------------------------------------------------------------------------------------------------------------------- 220 ! 221 222 !Program variables and types (finite elasticity part) 223 224 !Test program parameters 225 226 INTEGER(CMISSIntg) :: SolidMeshComponenetNumber 227 228 INTEGER(CMISSIntg), PARAMETER :: FieldGeometrySolidUserNumber=1 229 INTEGER(CMISSIntg), PARAMETER :: FieldGeometrySolidNumberOfVariables=1 230 INTEGER(CMISSIntg), PARAMETER :: FieldGeometrySolidNumberOfComponents=3 231 232 INTEGER(CMISSIntg), PARAMETER :: FieldFibreSolidUserNumber=2 233 INTEGER(CMISSIntg), PARAMETER :: FieldFibreSolidNumberOfVariables=1 234 INTEGER(CMISSIntg), PARAMETER :: FieldFibreSolidNumberOfComponents=3 235 236 INTEGER(CMISSIntg), PARAMETER :: FieldMaterialSolidUserNumber=3 237 INTEGER(CMISSIntg), PARAMETER :: FieldMaterialSolidNumberOfVariables=1 238 INTEGER(CMISSIntg), PARAMETER :: FieldMaterialSolidNumberOfComponents=3 239 240 INTEGER(CMISSIntg), PARAMETER :: FieldDependentSolidUserNumber=4 241 INTEGER(CMISSIntg), PARAMETER :: FieldDependentSolidNumberOfVariables=4 242 INTEGER(CMISSIntg), PARAMETER :: FieldDependentSolidNumberOfComponents=4 243 INTEGER(CMISSIntg), PARAMETER :: FieldDependentFluidNumberOfComponents=4 !(u,v,w,m) 244 245 INTEGER(CMISSIntg), PARAMETER :: EquationSetSolidUserNumber=1 246 INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldSolidUserNumber=25 247 248 INTEGER(CMISSIntg), PARAMETER :: SolidDisplMeshComponentNumber=1 249 INTEGER(CMISSIntg), PARAMETER :: SolidLagrMultMeshComponentNumber=2 250 INTEGER(CMISSIntg), PARAMETER :: SolidGeometryMeshComponentNumber=SolidDisplMeshComponentNumber 251 252 INTEGER(CMISSIntg), PARAMETER :: DarcyVelMeshComponentNumber=SolidLagrMultMeshComponentNumber 253 INTEGER(CMISSIntg), PARAMETER :: DarcyMassIncreaseMeshComponentNumber=SolidLagrMultMeshComponentNumber 254! INTEGER(CMISSIntg), PARAMETER :: DarcyGeometryMeshComponentNumber=SolidDisplMeshComponentNumber 255 256 INTEGER(CMISSIntg), PARAMETER :: GeneratedMeshUserNumber=32 257 !Program types 258 !Program variables 259 260 REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_START_TIME 261 REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_STOP_TIME 262 REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_THETA 263 REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_TIME_INCREMENT 264 265 !CMISS variables 266 267 TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsSolid 268 TYPE(CMISSEquationsType) :: EquationsSolid 269 TYPE(CMISSEquationsSetType) :: EquationsSetSolid 270 TYPE(CMISSFieldType) :: GeometricFieldSolid,FibreFieldSolid,MaterialFieldSolid 271 TYPE(CMISSFieldType) :: DependentFieldSolid,EquationsSetFieldSolid 272 TYPE(CMISSSolverType) :: SolverSolid 273 TYPE(CMISSSolverEquationsType) :: SolverEquationsSolid 274 275 !End - Program variables and types (finite elasticity part) 276 277 ! 278 !-------------------------------------------------------------------------------------------------------------------------------- 279 ! 280 281 282#ifdef WIN32 283 !Initialise QuickWin 284 QUICKWIN_WINDOW_CONFIG%TITLE="General Output" !Window title 285 QUICKWIN_WINDOW_CONFIG%NUMTEXTROWS=-1 !Max possible number of rows 286 QUICKWIN_WINDOW_CONFIG%MODE=QWIN$SCROLLDOWN 287 !Set the window parameters 288 QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG) 289 !If attempt fails set with system estimated values 290 IF(.NOT.QUICKWIN_STATUS) QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG) 291#endif 292 293 ! 294 !================================================================================================================================ 295 ! 296 NUMBER_GLOBAL_X_ELEMENTS=1 297 NUMBER_GLOBAL_Y_ELEMENTS=1 298 NUMBER_GLOBAL_Z_ELEMENTS=3 299 300 IF(NUMBER_GLOBAL_Z_ELEMENTS==0)THEN 301 NUMBER_OF_DIMENSIONS=2 302 ELSE 303 NUMBER_OF_DIMENSIONS=3 304 ENDIF 305 !PROBLEM CONTROL PANEL 306 307! BASIS_XI_INTERPOLATION_SOLID=CMISS_BASIS_LINEAR_LAGRANGE_INTERPOLATION 308 BASIS_XI_INTERPOLATION_SOLID=CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION 309 !Set geometric tolerance 310 GEOMETRY_TOLERANCE = 1.0E-12_CMISSDP 311 !Set initial values 312 INITIAL_FIELD_DARCY(1)=0.0_CMISSDP 313 INITIAL_FIELD_DARCY(2)=0.0_CMISSDP 314 INITIAL_FIELD_DARCY(3)=0.0_CMISSDP 315 INITIAL_FIELD_DARCY(4)=0.0_CMISSDP 316 !Set material parameters 317 POROSITY_PARAM_DARCY=0.1_CMISSDP 318 PERM_OVER_VIS_PARAM_DARCY=1.0_CMISSDP 319 !Set output parameter 320 !(NoOutput/ProgressOutput/TimingOutput/SolverOutput/SolverMatrixOutput) 321 DYNAMIC_SOLVER_DARCY_OUTPUT_TYPE=CMISS_SOLVER_PROGRESS_OUTPUT 322 LINEAR_SOLVER_DARCY_OUTPUT_TYPE=CMISS_SOLVER_SOLVER_OUTPUT 323 !(NoOutput/TimingOutput/MatrixOutput/ElementOutput) 324 EQUATIONS_DARCY_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT 325 326 !Set time parameter 327 DYNAMIC_SOLVER_DARCY_START_TIME=0.0_CMISSDP 328 DYNAMIC_SOLVER_DARCY_TIME_INCREMENT=1.0e-3_CMISSDP 329 DYNAMIC_SOLVER_DARCY_STOP_TIME=2_CMISSIntg * DYNAMIC_SOLVER_DARCY_TIME_INCREMENT 330 DYNAMIC_SOLVER_DARCY_THETA=1.0_CMISSDP !2.0_CMISSDP/3.0_CMISSDP 331 !Set result output parameter 332 DYNAMIC_SOLVER_DARCY_OUTPUT_FREQUENCY=1 333 !Set solver parameters 334 LINEAR_SOLVER_DARCY_DIRECT_FLAG=.TRUE. 335 RELATIVE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-05_CMISSDP 336 ABSOLUTE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-10_CMISSDP 337 DIVERGENCE_TOLERANCE=1.0E5_CMISSDP !default: 1.0E5 338 MAXIMUM_ITERATIONS=10000_CMISSIntg !default: 100000 339 RESTART_VALUE=30_CMISSIntg !default: 30 340 LINESEARCH_ALPHA=1.0_CMISSDP 341 342 343 ! 344 !================================================================================================================================ 345 ! 346 347 !INITIALISE OPENCMISS 348 349 CALL CMISSInitialise(WorldCoordinateSystem,WorldRegion,Err) 350 351 CALL CMISSErrorHandlingModeSet(CMISS_ERRORS_TRAP_ERROR,Err) 352 353 ! 354 !================================================================================================================================ 355 ! 356 357 !Set diagnostics 358 359 DIAG_LEVEL_LIST(1)=1 360 DIAG_LEVEL_LIST(2)=2 361 DIAG_LEVEL_LIST(3)=3 362 DIAG_LEVEL_LIST(4)=4 363 DIAG_LEVEL_LIST(5)=5 364 365 DIAG_ROUTINE_LIST(1)="WRITE_IP_INFO" 366! DIAG_ROUTINE_LIST(2)="FINITE_ELASTICITY_GAUSS_CAUCHY_TENSOR" 367 368 !CMISS_ALL_DIAG_TYPE/CMISS_IN_DIAG_TYPE/CMISS_FROM_DIAG_TYPE 369 CALL CMISSDiagnosticsSetOn(CMISS_IN_DIAG_TYPE,DIAG_LEVEL_LIST,"Diagnostics",DIAG_ROUTINE_LIST,Err) 370 371 !CMISS_ALL_TIMING_TYPE/CMISS_IN_TIMING_TYPE/CMISS_FROM_TIMING_TYPE 372 !TIMING_ROUTINE_LIST(1)="PROBLEM_FINITE_ELEMENT_CALCULATE" 373 !CALL TIMING_SET_ON(IN_TIMING_TYPE,.TRUE.,"",TIMING_ROUTINE_LIST,ERR,ERROR,*999) 374 375 ! 376 !================================================================================================================================ 377 ! 378 379 !Get the number of computational nodes and this computational node number 380 CALL CMISSComputationalNumberOfNodesGet(NumberOfComputationalNodes,Err) 381 CALL CMISSComputationalNodeNumberGet(ComputationalNodeNumber,Err) 382 383 NumberOfDomains = NumberOfComputationalNodes 384 write(*,*) "NumberOfDomains = ",NumberOfDomains 385 386 ! 387 !================================================================================================================================ 388 ! 389 390 !COORDINATE SYSTEM 391 392 CALL CMISSCoordinateSystem_Initialise(CoordinateSystem,Err) 393 CALL CMISSCoordinateSystem_CreateStart(CoordinateSystemUserNumber,CoordinateSystem,Err) 394 CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,NUMBER_OF_DIMENSIONS,Err) 395 CALL CMISSCoordinateSystem_CreateFinish(CoordinateSystem,Err) 396 397 ! 398 !================================================================================================================================ 399 ! 400 401 !REGION 402 !For a volume-coupled problem, solid and fluid are based in the same region 403 404 CALL CMISSRegion_Initialise(Region,Err) 405 CALL CMISSRegion_CreateStart(RegionUserNumber,WorldRegion,Region,Err) 406 CALL CMISSRegion_CoordinateSystemSet(Region,CoordinateSystem,Err) 407 CALL CMISSRegion_CreateFinish(Region,Err) 408 409 ! 410 !================================================================================================================================ 411 ! 412 413call READ_MESH('input/example-mesh-3els',MeshUserNumber,Region, Mesh,Bases,Nodes,Elements) 414! !BASES 415! !Define basis functions 416! CALL CMISSBasis_Initialise(LinearBasis,Err) 417! CALL CMISSBasis_CreateStart(LinearBasisUserNumber,LinearBasis,Err) 418! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(LinearBasis, & 419! & (/CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME/),Err) 420! !CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(LinearBasis,.true.,Err) 421! CALL CMISSBasis_CreateFinish(LinearBasis,Err) 422! 423! CALL CMISSBasis_Initialise(QuadraticBasis,Err) 424! CALL CMISSBasis_CreateStart(QuadraticBasisUserNumber,QuadraticBasis,Err) 425! CALL CMISSBasis_InterpolationXiSet(QuadraticBasis,(/CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION, & 426! & CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION,CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION/),Err) 427! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(QuadraticBasis, & 428! & (/CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME/),Err) 429! !CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(QuadraticBasis,.true.,Err) 430! CALL CMISSBasis_CreateFinish(QuadraticBasis,Err) 431! 432! CALL CMISSBasis_Initialise(CubicBasis,Err) 433! CALL CMISSBasis_CreateStart(CubicBasisUserNumber,CubicBasis,Err) 434! CALL CMISSBasis_InterpolationXiSet(CubicBasis,(/CMISS_BASIS_CUBIC_LAGRANGE_INTERPOLATION, & 435! & CMISS_BASIS_CUBIC_LAGRANGE_INTERPOLATION,CMISS_BASIS_CUBIC_LAGRANGE_INTERPOLATION/),Err) 436! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(CubicBasis, & 437! & (/CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME/),Err) 438! !CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(CubicBasis,.true.,Err) !Enable 3D interpolation on faces 439! CALL CMISSBasis_CreateFinish(CubicBasis,Err) 440! 441! !LinearBasis/QuadraticBasis/CubicBasis 442! Bases(1)=QuadraticBasis 443! Bases(2)=LinearBasis 444! 445! !Start the creation of a generated mesh in the region 446! CALL CMISSGeneratedMesh_Initialise(GeneratedMesh,Err) 447! CALL CMISSGeneratedMesh_CreateStart(GeneratedMeshUserNumber,Region,GeneratedMesh,Err) 448! CALL CMISSGeneratedMesh_TypeSet(GeneratedMesh,CMISS_GENERATED_MESH_REGULAR_MESH_TYPE,Err) 449! CALL CMISSGeneratedMesh_BasisSet(GeneratedMesh,Bases,Err) 450! IF(NUMBER_GLOBAL_Z_ELEMENTS==0) THEN 451! CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/X_DIM,Y_DIM/),Err) 452! CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS/),Err) 453! ELSE 454! CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/X_DIM,Y_DIM,Z_DIM/),Err) 455! CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS, & 456! & NUMBER_GLOBAL_Z_ELEMENTS/),Err) 457! ENDIF 458! CALL CMISSMesh_Initialise(Mesh,Err) 459! CALL CMISSGeneratedMesh_CreateFinish(GeneratedMesh,MeshUserNumber,Mesh,Err) 460! ------------------------------------------------------------------------------- 461 462 463 !GEOMETRIC FIELD 464 465 !Create a decomposition: 466 CALL CMISSDecomposition_Initialise(Decomposition,Err) 467 CALL CMISSDecomposition_CreateStart(DecompositionUserNumber,Mesh,Decomposition,Err) 468 !Set the decomposition to be a general decomposition with the specified number of domains 469 CALL CMISSDecomposition_TypeSet(Decomposition,CMISS_DECOMPOSITION_CALCULATED_TYPE,Err) 470 CALL CMISSDecomposition_NumberOfDomainsSet(Decomposition,NumberOfDomains,Err) 471 CALL CMISSDecomposition_CreateFinish(Decomposition,Err) 472 473 CALL CMISSField_Initialise(GeometricField,Err) 474 CALL CMISSField_CreateStart(GeometricFieldUserNumber,Region,GeometricField,Err) 475 CALL CMISSField_MeshDecompositionSet(GeometricField,Decomposition,Err) 476 CALL CMISSField_TypeSet(GeometricField,CMISS_FIELD_GEOMETRIC_TYPE,Err) 477 CALL CMISSField_NumberOfVariablesSet(GeometricField,FieldGeometryNumberOfVariables,Err) 478 CALL CMISSField_NumberOfComponentsSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,3,Err) 479 CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidGeometryMeshComponentNumber,Err) 480 CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidGeometryMeshComponentNumber,Err) 481 CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidGeometryMeshComponentNumber,Err) 482 CALL CMISSField_CreateFinish(GeometricField,Err) 483CALL READ_NODES('input/example-nodes-3els',GeometricField) 484! CALL CMISSGeneratedMesh_GeometricParametersCalculate(GeneratedMesh,GeometricField,Err) 485 486 !-------------------------------------------------------------------------------------------------------------------------------- 487 ! Solid 488 489 !Create a decomposition 490 491 !Create a field to put the geometry (defualt is geometry) 492 493 SolidMeshComponenetNumber = SolidGeometryMeshComponentNumber 494 495 CALL CMISSField_Initialise(GeometricFieldSolid,Err) 496 CALL CMISSField_CreateStart(FieldGeometrySolidUserNumber,Region,GeometricFieldSolid,Err) 497 CALL CMISSField_MeshDecompositionSet(GeometricFieldSolid,Decomposition,Err) 498 CALL CMISSField_TypeSet(GeometricFieldSolid,CMISS_FIELD_GEOMETRIC_TYPE,Err) 499 CALL CMISSField_NumberOfVariablesSet(GeometricFieldSolid,FieldGeometrySolidNumberOfVariables,Err) 500 CALL CMISSField_NumberOfComponentsSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldGeometrySolidNumberOfComponents,Err) 501 CALL CMISSField_ComponentMeshComponentSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidMeshComponenetNumber,Err) 502 CALL CMISSField_ComponentMeshComponentSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidMeshComponenetNumber,Err) 503 CALL CMISSField_ComponentMeshComponentSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidMeshComponenetNumber,Err) 504 CALL CMISSField_CreateFinish(GeometricFieldSolid,Err) 505 !Set the mesh component to be used by the field components. 506CALL READ_NODES('input/example-nodes-3els',GeometricFieldSolid) 507! CALL CMISSGeneratedMesh_GeometricParametersCalculate(GeneratedMesh,GeometricFieldSolid,Err) 508 509 !Create a fibre field and attach it to the geometric field 510 CALL CMISSField_Initialise(FibreFieldSolid,Err) 511 CALL CMISSField_CreateStart(FieldFibreSolidUserNumber,Region,FibreFieldSolid,Err) 512 CALL CMISSField_TypeSet(FibreFieldSolid,CMISS_FIELD_FIBRE_TYPE,Err) 513 CALL CMISSField_MeshDecompositionSet(FibreFieldSolid,Decomposition,Err) 514 CALL CMISSField_GeometricFieldSet(FibreFieldSolid,GeometricFieldSolid,Err) 515 CALL CMISSField_NumberOfVariablesSet(FibreFieldSolid,FieldFibreSolidNumberOfVariables,Err) 516 CALL CMISSField_NumberOfComponentsSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldFibreSolidNumberOfComponents,Err) 517 CALL CMISSField_ComponentMeshComponentSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidGeometryMeshComponentNumber,Err) 518 CALL CMISSField_ComponentMeshComponentSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidGeometryMeshComponentNumber,Err) 519 CALL CMISSField_ComponentMeshComponentSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidGeometryMeshComponentNumber,Err) 520 CALL CMISSField_CreateFinish(FibreFieldSolid,Err) 521 522 ! end Solid 523 !-------------------------------------------------------------------------------------------------------------------------------- 524 525 ! 526 !================================================================================================================================ 527 ! 528 529 !EQUATIONS SETS 530 531 !Create the equations set for ALE Darcy 532 CALL CMISSField_Initialise(EquationsSetFieldDarcy,Err) 533 CALL CMISSEquationsSet_Initialise(EquationsSetDarcy,Err) 534 CALL CMISSEquationsSet_CreateStart(EquationsSetUserNumberDarcy,Region,GeometricField,CMISS_EQUATIONS_SET_FLUID_MECHANICS_CLASS, & 535 & CMISS_EQUATIONS_SET_DARCY_EQUATION_TYPE,CMISS_EQUATIONS_SET_INCOMPRESS_ELASTICITY_DRIVEN_DARCY_SUBTYPE,& 536 & EquationsSetFieldUserNumberDarcy,EquationsSetFieldDarcy,EquationsSetDarcy,Err) 537 CALL CMISSEquationsSet_CreateFinish(EquationsSetDarcy,Err) 538 539 !Create the equations set for the solid 540 CALL CMISSField_Initialise(EquationsSetFieldSolid,Err) 541 CALL CMISSEquationsSet_Initialise(EquationsSetSolid,Err) 542 CALL CMISSEquationsSet_CreateStart(EquationSetSolidUserNumber,Region,FibreFieldSolid,CMISS_EQUATIONS_SET_ELASTICITY_CLASS, & 543 & CMISS_EQUATIONS_SET_FINITE_ELASTICITY_TYPE,CMISS_EQUATIONS_SET_INCOMPRESS_ELASTICITY_DRIVEN_DARCY_SUBTYPE,& 544 & EquationsSetFieldSolidUserNumber,EquationsSetFieldSolid,EquationsSetSolid,Err) 545 CALL CMISSEquationsSet_CreateFinish(EquationsSetSolid,Err) 546 547 !-------------------------------------------------------------------------------------------------------------------------------- 548 ! Solid Materials Field 549 550 !Create a material field and attach it to the geometric field 551 CALL CMISSField_Initialise(MaterialFieldSolid,Err) 552 ! 553 CALL CMISSField_CreateStart(FieldMaterialSolidUserNumber,Region,MaterialFieldSolid,Err) 554 ! 555 CALL CMISSField_TypeSet(MaterialFieldSolid,CMISS_FIELD_MATERIAL_TYPE,Err) 556 CALL CMISSField_MeshDecompositionSet(MaterialFieldSolid,Decomposition,Err) 557 CALL CMISSField_GeometricFieldSet(MaterialFieldSolid,GeometricFieldSolid,Err) 558 CALL CMISSField_NumberOfVariablesSet(MaterialFieldSolid,FieldMaterialSolidNumberOfVariables,Err) 559 CALL CMISSField_NumberOfComponentsSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldMaterialSolidNumberOfComponents,Err) 560 CALL CMISSField_ComponentMeshComponentSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidGeometryMeshComponentNumber,Err) 561 CALL CMISSField_ComponentMeshComponentSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidGeometryMeshComponentNumber,Err) 562 CALL CMISSField_ComponentMeshComponentSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidGeometryMeshComponentNumber,Err) 563 ! 564 CALL CMISSField_CreateFinish(MaterialFieldSolid,Err) 565 566 !Set material parameters 567 CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1, & 568 & 2.0_CMISSDP,Err) 569! CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,2.0e3_CMISSDP,Err) 570 CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,2, & 571 & 6.0_CMISSDP,Err) 572! CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,2,33.0_CMISSDP,Err) 573 CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,3, & 574 & 10.0_CMISSDP,Err) 575 576 577 CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetSolid,FieldMaterialSolidUserNumber,MaterialFieldSolid,Err) 578 CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetSolid,Err) 579 580 ! end Solid 581 !-------------------------------------------------------------------------------------------------------------------------------- 582 583 584 ! 585 !================================================================================================================================ 586 ! 587 588 !DEPENDENT FIELDS 589 590 !-------------------------------------------------------------------------------------------------------------------------------- 591 ! Solid 592 593 !Create a dependent field with four variables (U, DelUDelN = solid, V, DelVDelN = Darcy) and four components 594 CALL CMISSField_Initialise(DependentFieldSolid,Err) 595 ! 596 CALL CMISSField_CreateStart(FieldDependentSolidUserNumber,Region,DependentFieldSolid,Err) 597 ! 598 CALL CMISSField_TypeSet(DependentFieldSolid,CMISS_FIELD_GENERAL_TYPE,Err) 599 CALL CMISSField_MeshDecompositionSet(DependentFieldSolid,Decomposition,Err) 600 CALL CMISSField_GeometricFieldSet(DependentFieldSolid,GeometricFieldSolid,Err) 601 CALL CMISSField_DependentTypeSet(DependentFieldSolid,CMISS_FIELD_DEPENDENT_TYPE,Err) 602 CALL CMISSField_NumberOfVariablesSet(DependentFieldSolid,FieldDependentSolidNumberOfVariables,Err) 603 CALL CMISSField_VariableTypesSet(DependentFieldSolid,(/CMISS_FIELD_U_VARIABLE_TYPE, & 604 & CMISS_FIELD_DELUDELN_VARIABLE_TYPE,CMISS_FIELD_V_VARIABLE_TYPE,CMISS_FIELD_DELVDELN_VARIABLE_TYPE/),Err) 605 CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldDependentSolidNumberOfComponents,Err) 606 CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE, & 607 & FieldDependentSolidNumberOfComponents,Err) 608 CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,FieldDependentFluidNumberOfComponents,Err) 609 CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE, & 610 & FieldDependentFluidNumberOfComponents,Err) 611 ! 612 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidDisplMeshComponentNumber,Err) 613 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidDisplMeshComponentNumber,Err) 614 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidDisplMeshComponentNumber,Err) 615 CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4, & 616 & CMISS_FIELD_NODE_BASED_INTERPOLATION, & 617 & Err) 618! CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4,CMISS_FIELD_ELEMENT_BASED_INTERPOLATION,Err) 619! CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4,SolidMeshComponenetNumber,Err) 620 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4,SolidLagrMultMeshComponentNumber,Err) 621 ! 622 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,1, & 623 & SolidDisplMeshComponentNumber, & 624 & Err) 625 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,2, & 626 & SolidDisplMeshComponentNumber, & 627 & Err) 628 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,3, & 629 & SolidDisplMeshComponentNumber, & 630 & Err) 631 CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4, & 632 & CMISS_FIELD_NODE_BASED_INTERPOLATION,Err) 633! CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4, & 634! & CMISS_FIELD_ELEMENT_BASED_INTERPOLATION,Err) 635! CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4,SolidMeshComponenetNumber,Err) 636 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4, & 637 & SolidLagrMultMeshComponentNumber, & 638 & Err) 639 640 !For this equation type, MESH_COMPONENT_NUMBER_PRESSURE is actually the mass increase component as the pressure is taken from the solid equations 641 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,DarcyVelMeshComponentNumber,Err) 642 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,2,DarcyVelMeshComponentNumber,Err) 643 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,3,DarcyVelMeshComponentNumber,Err) 644! CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,4,MESH_COMPONENT_NUMBER_PRESSURE,Err) 645 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,4, & 646 & DarcyMassIncreaseMeshComponentNumber, & 647 & Err) 648 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,1,DarcyVelMeshComponentNumber, & 649 & Err) 650 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,2,DarcyVelMeshComponentNumber, & 651 & Err) 652 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,3,DarcyVelMeshComponentNumber, & 653 & Err) 654! CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,4,MESH_COMPONENT_NUMBER_PRESSURE,Err) 655 CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,4, & 656 & DarcyMassIncreaseMeshComponentNumber,Err) 657 658 CALL CMISSField_CreateFinish(DependentFieldSolid,Err) 659 ! 660 CALL CMISSEquationsSet_DependentCreateStart(EquationsSetSolid,FieldDependentSolidUserNumber,DependentFieldSolid,Err) 661 CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetSolid,Err) 662 663 ! end Solid 664 !-------------------------------------------------------------------------------------------------------------------------------- 665 666 667 !Create the equations set dependent field variables for ALE Darcy 668 CALL CMISSEquationsSet_DependentCreateStart(EquationsSetDarcy,FieldDependentSolidUserNumber, & ! ??? UserNumber ??? 669 & DependentFieldSolid,Err) 670 CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetDarcy,Err) 671 672 !Initialise dependent field (velocity components,mass increase) 673 DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS+1 674 CALL CMISSField_ComponentValuesInitialise(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 675 & COMPONENT_NUMBER,INITIAL_FIELD_DARCY(COMPONENT_NUMBER),Err) 676 ENDDO 677 678 679 ! 680 !================================================================================================================================ 681 ! 682 683call READ_FIELD('input/example-field-darcy',MaterialsFieldUserNumberDarcy,Region,GeometricField,MaterialsFieldDarcy) 684CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetDarcy,MaterialsFieldUserNumberDarcy,MaterialsFieldDarcy,Err) 685CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetDarcy,Err) 686! !MATERIALS FIELDS 687! 688! !Create the equations set materials field variables for ALE Darcy 689! CALL CMISSField_Initialise(MaterialsFieldDarcy,Err) 690! CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetDarcy,MaterialsFieldUserNumberDarcy, & 691! & MaterialsFieldDarcy,Err) 692! !Finish the equations set materials field variables 693! CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetDarcy,Err) 694! CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDarcy,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 695! & MaterialsFieldUserNumberDarcyPorosity,POROSITY_PARAM_DARCY,Err) 696! CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDarcy,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 697! & MaterialsFieldUserNumberDarcyPermOverVis,PERM_OVER_VIS_PARAM_DARCY,Err) 698 699 ! 700 !================================================================================================================================ 701 ! 702 703 !EQUATIONS SET EQUATIONS 704 705 !Darcy 706 CALL CMISSEquations_Initialise(EquationsDarcy,Err) 707 CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetDarcy,EquationsDarcy,Err) 708 CALL CMISSEquations_SparsityTypeSet(EquationsDarcy,CMISS_EQUATIONS_SPARSE_MATRICES,Err) 709 CALL CMISSEquations_OutputTypeSet(EquationsDarcy,EQUATIONS_DARCY_OUTPUT,Err) 710 CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetDarcy,Err) 711 712 !Solid 713 CALL CMISSEquations_Initialise(EquationsSolid,Err) 714 CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetSolid,EquationsSolid,Err) 715 CALL CMISSEquations_SparsityTypeSet(EquationsSolid,CMISS_EQUATIONS_SPARSE_MATRICES,Err) 716 CALL CMISSEquations_OutputTypeSet(EquationsSolid,CMISS_EQUATIONS_NO_OUTPUT,Err) 717 CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetSolid,Err) 718 719 ! 720 !================================================================================================================================ 721 ! 722 723 !-------------------------------------------------------------------------------------------------------------------------------- 724 ! Solid 725 726 !Initialise dependent field from undeformed geometry and displacement bcs and set hydrostatic pressure 727 CALL CMISSField_ParametersToFieldParametersComponentCopy(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE, & 728 & CMISS_FIELD_VALUES_SET_TYPE, & 729 & 1,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,Err) 730 CALL CMISSField_ParametersToFieldParametersComponentCopy(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE, & 731 & CMISS_FIELD_VALUES_SET_TYPE, & 732 & 2,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,2,Err) 733 CALL CMISSField_ParametersToFieldParametersComponentCopy(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE, & 734 & CMISS_FIELD_VALUES_SET_TYPE, & 735 & 3,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,3,Err) 736 CALL CMISSField_ComponentValuesInitialise(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,4, & 737 & 0.0_CMISSDP, & 738 & Err) 739 740 741 ! 742 !================================================================================================================================ 743 ! 744 745 !PROBLEMS 746 747 CALL CMISSProblem_Initialise(Problem,Err) 748 CALL CMISSControlLoop_Initialise(ControlLoop,Err) 749 CALL CMISSProblem_CreateStart(ProblemUserNumber,Problem,Err) 750 CALL CMISSProblem_SpecificationSet(Problem,CMISS_PROBLEM_MULTI_PHYSICS_CLASS,CMISS_PROBLEM_FINITE_ELASTICITY_DARCY_TYPE, & 751 & CMISS_PROBLEM_QUASISTATIC_ELASTICITY_TRANSIENT_DARCY_SUBTYPE,Err) 752 CALL CMISSProblem_CreateFinish(Problem,Err) 753 754 CALL CMISSProblem_ControlLoopCreateStart(Problem,Err) 755 CALL CMISSProblem_ControlLoopGet(Problem,CMISS_CONTROL_LOOP_NODE,ControlLoop,Err) 756! CALL CMISSControlLoop_MaximumIterationsSet(ControlLoop,2,Err) 757 CALL CMISSControlLoop_TimesSet(ControlLoop,DYNAMIC_SOLVER_DARCY_START_TIME,DYNAMIC_SOLVER_DARCY_STOP_TIME, & 758 & DYNAMIC_SOLVER_DARCY_TIME_INCREMENT,Err) 759 CALL CMISSControlLoop_TimeOutputSet(ControlLoop,DYNAMIC_SOLVER_DARCY_OUTPUT_FREQUENCY,Err) 760! CALL CMISSControlLoop_OutputTypeSet(ControlLoop,CMISS_CONTROL_LOOP_PROGRESS_OUTPUT,Err) 761 CALL CMISSProblem_ControlLoopCreateFinish(Problem,Err) 762 763 ! 764 !================================================================================================================================ 765 ! 766 767 !SOLVERS 768 769 CALL CMISSSolver_Initialise(SolverSolid,Err) 770 CALL CMISSSolver_Initialise(DynamicSolverDarcy,Err) 771 CALL CMISSSolver_Initialise(LinearSolverDarcy,Err) 772 773 CALL CMISSProblem_SolversCreateStart(Problem,Err) 774 775 ! Solid 776 CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopSolidNumber,CMISS_CONTROL_LOOP_NODE/), & 777 & SolverSolidIndex,SolverSolid,Err) 778 CALL CMISSSolver_OutputTypeSet(SolverSolid,CMISS_SOLVER_PROGRESS_OUTPUT,Err) 779! CALL CMISSSolver_NewtonJacobianCalculationTypeSet(SolverSolid,CMISS_SOLVER_NEWTON_JACOBIAN_FD_CALCULATED,Err) 780 CALL CMISSSolver_NewtonJacobianCalculationTypeSet(SolverSolid,CMISS_SOLVER_NEWTON_JACOBIAN_EQUATIONS_CALCULATED,Err) 781 782 CALL CMISSSolver_NewtonAbsoluteToleranceSet(SolverSolid,ABSOLUTE_TOLERANCE,Err) 783 CALL CMISSSolver_NewtonRelativeToleranceSet(SolverSolid,RELATIVE_TOLERANCE,Err) 784 CALL CMISSSolver_NewtonMaximumIterationsSet(SolverSolid,MAXIMUM_ITERATIONS,Err) 785 786! CALL CMISSSolverNonLinearTypeSet(SolverSolid,CMISS_SOLVER_NONLINEAR_NEWTON,Err) 787! CALL CMISSSolver_LibraryTypeSet(SolverSolid,CMISS_SOLVER_PETSC_LIBRARY,Err) 788 789! CALL CMISSSolver_NewtonLinearSolverGet(SolverSolid,LinearSolverSolid,Err) 790! CALL CMISSSolver_LinearTypeSet(LinearSolverSolid,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err) 791 792 793 !Darcy 794 CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopFluidNumber,CMISS_CONTROL_LOOP_NODE/), & 795 & SolverDarcyIndex,DynamicSolverDarcy,Err) 796 CALL CMISSSolver_OutputTypeSet(DynamicSolverDarcy,DYNAMIC_SOLVER_DARCY_OUTPUT_TYPE,Err) 797 CALL CMISSSolver_DynamicThetaSet(DynamicSolverDarcy,DYNAMIC_SOLVER_DARCY_THETA,Err) 798! CALL CMISSSolverDynamicDynamicSet(DynamicSolverDarcy,.TRUE.,Err) 799 CALL CMISSSolver_DynamicLinearSolverGet(DynamicSolverDarcy,LinearSolverDarcy,Err) 800 IF(LINEAR_SOLVER_DARCY_DIRECT_FLAG) THEN 801 CALL CMISSSolver_LinearTypeSet(LinearSolverDarcy,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err) 802 CALL CMISSSolver_LibraryTypeSet(LinearSolverDarcy,CMISS_SOLVER_MUMPS_LIBRARY,Err) 803 ELSE 804 CALL CMISSSolver_LinearTypeSet(LinearSolverDarcy,CMISS_SOLVER_LINEAR_ITERATIVE_SOLVE_TYPE,Err) 805 CALL CMISSSolver_LinearIterativeMaximumIterationsSet(LinearSolverDarcy,MAXIMUM_ITERATIONS,Err) 806 CALL CMISSSolver_LinearIterativeDivergenceToleranceSet(LinearSolverDarcy,DIVERGENCE_TOLERANCE,Err) 807 CALL CMISSSolver_LinearIterativeRelativeToleranceSet(LinearSolverDarcy,RELATIVE_TOLERANCE,Err) 808 CALL CMISSSolver_LinearIterativeAbsoluteToleranceSet(LinearSolverDarcy,ABSOLUTE_TOLERANCE,Err) 809 CALL CMISSSolver_LinearIterativeGMRESRestartSet(LinearSolverDarcy,RESTART_VALUE,Err) 810 ENDIF 811 812 CALL CMISSProblem_SolversCreateFinish(Problem,Err) 813 814 ! 815 !================================================================================================================================ 816 ! 817 818 !SOLVER EQUATIONS 819 820 CALL CMISSSolver_Initialise(SolverSolid,Err) 821 CALL CMISSSolver_Initialise(LinearSolverDarcy,Err) 822 823 CALL CMISSSolverEquations_Initialise(SolverEquationsSolid,Err) 824 CALL CMISSSolverEquations_Initialise(SolverEquationsDarcy,Err) 825 826 CALL CMISSProblem_SolverEquationsCreateStart(Problem,Err) 827 ! 828 !Get the finite elasticity solver equations 829 CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopSolidNumber,CMISS_CONTROL_LOOP_NODE/), & 830 & SolverSolidIndex,SolverSolid,Err) 831 CALL CMISSSolver_SolverEquationsGet(SolverSolid,SolverEquationsSolid,Err) 832 CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsSolid,CMISS_SOLVER_SPARSE_MATRICES,Err) 833 CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsSolid,EquationsSetSolid,EquationsSetIndex,Err) 834 ! 835 !Get the Darcy solver equations 836 CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopFluidNumber,CMISS_CONTROL_LOOP_NODE/), & 837 & SolverDarcyIndex,LinearSolverDarcy,Err) 838 CALL CMISSSolver_SolverEquationsGet(LinearSolverDarcy,SolverEquationsDarcy,Err) 839 CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsDarcy,CMISS_SOLVER_SPARSE_MATRICES,Err) 840 CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsDarcy,EquationsSetDarcy,EquationsSetIndex,Err) 841 ! 842 CALL CMISSProblem_SolverEquationsCreateFinish(Problem,Err) 843 844 ! 845 !================================================================================================================================ 846 ! 847 848 !------------------------------------ 849 ! ASSIGN BOUNDARY CONDITIONS - SOLID (absolute nodal parameters) 850 !Solid is computed in absolute position, rather than displacement. Thus BCs for absolute position 851 CALL CMISSBoundaryConditions_Initialise(BoundaryConditionsSolid,Err) 852 CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsSolid,BoundaryConditionsSolid,Err) 853 854! !Get surfaces 855! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_FRONT_SURFACE, & 856! & Face1Nodes,FaceXi(1),Err) 857! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BACK_SURFACE, & 858! & Face2Nodes,FaceXi(2),Err) 859! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_RIGHT_SURFACE, & 860! & Face3Nodes,FaceXi(3),Err) 861! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_LEFT_SURFACE, & 862! & Face4Nodes,FaceXi(4),Err) 863! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_TOP_SURFACE, & 864! & Face5Nodes,FaceXi(5),Err) 865! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BOTTOM_SURFACE, & 866! & Face6Nodes,FaceXi(6),Err) 867 868! ! write(*,'("Face1Nodes=[",100I3,"]")') Face1Nodes 869! ! write(*,'("Face2Nodes=[",100I3,"]")') Face2Nodes 870! ! write(*,'("Face3Nodes=[",100I3,"]")') Face3Nodes 871! ! write(*,'("Face4Nodes=[",100I3,"]")') Face4Nodes 872! ! write(*,'("Face5Nodes=[",100I3,"]")') Face5Nodes 873! ! write(*,'("Face6Nodes=[",100I3,"]")') Face6Nodes 874! ! write(*,'("FaceXi=[",100I3,"]")') FaceXi 875! 876allocate(Face1Nodes(21),Face2Nodes(21),Face3Nodes(21),Face4Nodes(21),Face5Nodes(9),Face6Nodes(9)) 877Face1Nodes=[ 1, 17, 2, 22, 23, 24, 5, 31, 6, 36, 37, 38, 9, 45, 10, 50, 51, 52, 13, 59, 14] 878Face2Nodes=[ 3, 21, 4, 28, 29, 30, 7, 35, 8, 42, 43, 44, 11, 49, 12, 56, 57, 58, 15, 63, 16] 879Face3Nodes=[ 2, 20, 4, 24, 27, 30, 6, 34, 8, 38, 41, 44, 10, 48, 12, 52, 55, 58, 14, 62, 16] 880Face4Nodes=[ 1, 18, 3, 22, 25, 28, 5, 32, 7, 36, 39, 42, 9, 46, 11, 50, 53, 56, 13, 60, 15] 881Face5Nodes=[ 13, 59, 14, 60, 61, 62, 15, 63, 16] 882Face6Nodes=[ 1, 17, 2, 18, 19, 20, 3, 21, 4] 883FaceXi=[ -2, 2, 1, -1, 3, -3] 884 885 886 ! Fix the bottom in z direction 887 DO NN=1,SIZE(Face6Nodes,1) 888 NODE=Face6Nodes(NN) 889 CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err) 890 IF(NodeDomain==ComputationalNodeNumber) THEN 891 CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,3, & 892 & ZCoord,Err) 893 CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,3, & 894 & CMISS_BOUNDARY_CONDITION_FIXED,ZCoord,Err) 895 WRITE(*,*) "FIXING NODE",NODE,"AT BOTTOM IN Z DIRECTION" 896 ENDIF 897 ENDDO 898 899 ! Fix the top in z direction 900 DO NN=1,SIZE(Face5Nodes,1) 901 NODE=Face5Nodes(NN) 902 CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err) 903 IF(NodeDomain==ComputationalNodeNumber) THEN 904 CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,3, & 905 & ZCoord,Err) 906 CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,3, & 907 & CMISS_BOUNDARY_CONDITION_FIXED,ZCoord,Err) 908 WRITE(*,*) "FIXING NODE",NODE,"AT TOP IN Z DIRECTION" 909 ENDIF 910 ENDDO 911 912 !Fix more nodes at the bottom to stop free body motion 913 X_FIXED=.FALSE. 914 Y_FIXED=.FALSE. 915 DO NN=1,SIZE(Face6Nodes,1) 916 NODE=Face6Nodes(NN) 917 CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err) 918 IF(NodeDomain==ComputationalNodeNumber) THEN 919 CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,1, & 920 & XCoord,Err) 921 CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,2, & 922 & YCoord,Err) 923 924 !Fix Origin displacement in x and y (z already fixed) 925 IF(ABS(XCoord)<1.0E-6_CMISSDP.AND.ABS(YCoord)<1.0E-6_CMISSDP) THEN 926 CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,1, & 927 & CMISS_BOUNDARY_CONDITION_FIXED,XCoord,Err) 928 CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,2, & 929 & CMISS_BOUNDARY_CONDITION_FIXED,YCoord,Err) 930 WRITE(*,*) "FIXING ORIGIN NODE",NODE,"IN X AND Y DIRECTION" 931 X_FIXED=.TRUE. 932 Y_FIXED=.TRUE. 933 ENDIF 934 935 !Fix nodal displacements at (X_DIM,0) in y 936 IF(ABS(XCoord - X_DIM)<1.0E-6_CMISSDP .AND. ABS(YCoord)<1.0E-6_CMISSDP) THEN 937 CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,2, & 938 & CMISS_BOUNDARY_CONDITION_FIXED,YCoord,Err) 939 WRITE(*,*) "FIXING NODES",NODE,"AT (X_DIM,0) IN Y DIRECTION" 940 Y_FIXED=.TRUE. 941 ENDIF 942 943 !Fix nodal displacements at (0,Y_DIM) in x 944 IF(ABS(XCoord)<1.0E-6_CMISSDP .AND. ABS(YCoord - Y_DIM)<1.0E-6_CMISSDP) THEN 945 CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,1, & 946 & CMISS_BOUNDARY_CONDITION_FIXED,XCoord,Err) 947 WRITE(*,*) "FIXING NODES",NODE,"AT (0,Y_DIM) IN X DIRECTION" 948 X_FIXED=.TRUE. 949 ENDIF 950 951 ENDIF 952 ENDDO 953! CALL MPI_REDUCE(X_FIXED,X_OKAY,1,MPI_LOGICAL,MPI_LOR,0,MPI_COMM_WORLD,MPI_IERROR) 954! CALL MPI_REDUCE(Y_FIXED,Y_OKAY,1,MPI_LOGICAL,MPI_LOR,0,MPI_COMM_WORLD,MPI_IERROR) 955! IF(ComputationalNodeNumber==0) THEN 956! IF(.NOT.(X_OKAY.AND.Y_OKAY)) THEN 957! WRITE(*,*) "Free body motion could not be prevented!" 958! CALL CMISSFinalise(Err) 959! STOP 960! ENDIF 961! ENDIF 962 963 CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquationsSolid,Err) 964 !------------------------------------ 965 966 967 !------------------------------------ 968 ! ASSIGN BOUNDARY CONDITIONS - FLUID 969 CALL CMISSBoundaryConditions_Initialise(BoundaryConditionsDarcy,Err) 970 CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsDarcy,BoundaryConditionsDarcy,Err) 971 972! !Get surfaces 973! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_FRONT_SURFACE, & 974! & Face7Nodes,FaceXi(1),Err) 975! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BACK_SURFACE, & 976! & Face8Nodes,FaceXi(2),Err) 977! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_RIGHT_SURFACE, & 978! & Face9Nodes,FaceXi(3),Err) 979! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_LEFT_SURFACE, & 980! & Face10Nodes,FaceXi(4),Err) 981! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_TOP_SURFACE, & 982! & Face11Nodes,FaceXi(5),Err) 983! CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BOTTOM_SURFACE, & 984! & Face12Nodes,FaceXi(6),Err) 985 986! ! write(*,'("Face7Nodes=[",100I3,"]")') Face7Nodes 987! ! write(*,'("Face8Nodes=[",100I3,"]")') Face8Nodes 988! ! write(*,'("Face9Nodes=[",100I3,"]")') Face9Nodes 989!…
Large files files are truncated, but you can click here to view the full file