PageRenderTime 66ms CodeModel.GetById 21ms RepoModel.GetById 0ms app.codeStats 0ms

/MultiPhysics/Poroelasticity/FiniteElasticityDarcy/IncompressibleElasticityDrivenDarcyIO/src/IncompressibleElasticityDrivenDarcyIOExample.f90

http://github.com/xyan075/examples
FORTRAN Modern | 1137 lines | 571 code | 184 blank | 382 comment | 0 complexity | fa0189050314fecf43370947ff8449ab MD5 | raw file
    

  1. !> \file
    2. 

  2. !> \author Christian Michler, Adam Reeve
    3. 

  3. !> \brief This is an example program to solve a coupled Finite Elastiticity Darcy equation using openCMISS calls.
    4. 

  4. !>
    5. 

  5. !> \section LICENSE
    6. 

  6. !>
    7. 

  7. !> Version: MPL 1.1/GPL 2.0/LGPL 2.1
    8. 

  8. !>
    9. 

  9. !> The contents of this file are subject to the Mozilla Public License
    10. 

  10. !> Version 1.1 (the "License"); you may not use this file except in
    11. 

  11. !> compliance with the License. You may obtain a copy of the License at
    12. 

  12. !> http://www.mozilla.org/MPL/
    13. 

  13. !>
    14. 

  14. !> Software distributed under the License is distributed on an "AS IS"
    15. 

  15. !> basis, WITHOUT WARRANTY OF ANY KIND, either express or implied. See the
    16. 

  16. !> License for the specific language governing rights and limitations
    17. 

  17. !> under the License.
    18. 

  18. !>
    19. 

  19. !> The Original Code is OpenCMISS
    20. 

  20. !>
    21. 

  21. !> The Initial Developer of the Original Code is University of Auckland,
    22. 

  22. !> Auckland, New Zealand and University of Oxford, Oxford, United
    23. 

  23. !> Kingdom. Portions created by the University of Auckland and University
    24. 

  24. !> of Oxford are Copyright (C) 2007 by the University of Auckland and
    25. 

  25. !> the University of Oxford. All Rights Reserved.
    26. 

  26. !>
    27. 

  27. !> Contributor(s):
    28. 

  28. !>
    29. 

  29. !> Alternatively, the contents of this file may be used under the terms of
    30. 

  30. !> either the GNU General Public License Version 2 or later (the "GPL"), or
    31. 

  31. !> the GNU Lesser General Public License Version 2.1 or later (the "LGPL"),
    32. 

  32. !> in which case the provisions of the GPL or the LGPL are applicable instead
    33. 

  33. !> of those above. If you wish to allow use of your version of this file only
    34. 

  34. !> under the terms of either the GPL or the LGPL, and not to allow others to
    35. 

  35. !> use your version of this file under the terms of the MPL, indicate your
    36. 

  36. !> decision by deleting the provisions above and replace them with the notice
    37. 

  37. !> and other provisions required by the GPL or the LGPL. If you do not delete
    38. 

  38. !> the provisions above, a recipient may use your version of this file under
    39. 

  39. !> the terms of any one of the MPL, the GPL or the LGPL.
    40. 

  40. !>
    41. 

  41. 

  42. !> \example MultiPhysics/Poroelasticity/FiniteElasticityDarcy/IncompressibleElasticityDrivenDarcy/src/IncompressibleElasticityDrivenDarcyExample.f90
    43. 

  43. !! Example program to solve coupled FiniteElasticityDarcy equations using OpenCMISS calls.
    44. 

  44. !! \par Latest Builds:
    45. 

  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. 

  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. 

  47. !!
    48. 

  48. !<
    49. 

  49. 

  50. ! !
    51. 

  51. ! !  This example considers a coupled Finite Elasticity Darcy problem
    52. 

  52. ! !
    53. 

  53. 

  54. !> Main program
    55. 

  55. 

  56. PROGRAM FINITEELASTICITYDARCYIOEXAMPLE
    57. 

  57. 

  58.   !
    59. 

  59.   !================================================================================================================================
    60. 

  60.   !
    61. 

  61. 

  62.   !PROGRAM LIBRARIES
    63. 

  63. 

  64.   USE OPENCMISS
    65. 

  65. !   USE FLUID_MECHANICS_IO_ROUTINES
    66. 

  66.   USE IOSTUFF
    67. 

  67.   USE MPI
    68. 

  68. 

  69. #ifdef WIN32
    70. 

  70.   USE IFQWINCMISS
    71. 

  71. #endif
    72. 

  72. 

  73.   !
    74. 

  74.   !================================================================================================================================
    75. 

  75.   !
    76. 

  76. 

  77.   !PROGRAM VARIABLES AND TYPES
    78. 

  78. 

  79.   IMPLICIT NONE
    80. 

  80. 

  81.   !Test program parameters
    82. 

  82. 

  83.   REAL(CMISSDP), PARAMETER :: Y_DIM=1.0_CMISSDP
    84. 

  84.   REAL(CMISSDP), PARAMETER :: X_DIM=1.0_CMISSDP
    85. 

  85.   REAL(CMISSDP), PARAMETER :: Z_DIM=1.0_CMISSDP
    86. 

  86. 

  87.   INTEGER(CMISSIntg), PARAMETER :: LinearBasisUserNumber=1
    88. 

  88.   INTEGER(CMISSIntg), PARAMETER :: QuadraticBasisUserNumber=2
    89. 

  89.   INTEGER(CMISSIntg), PARAMETER :: CubicBasisUserNumber=3
    90. 

  90. 

  91.   INTEGER(CMISSIntg), PARAMETER :: CoordinateSystemUserNumber=1
    92. 

  92.   INTEGER(CMISSIntg), PARAMETER :: RegionUserNumber=2
    93. 

  93.   INTEGER(CMISSIntg), PARAMETER :: MeshUserNumber=3
    94. 

  94.   INTEGER(CMISSIntg), PARAMETER :: DecompositionUserNumber=4
    95. 

  95.   INTEGER(CMISSIntg), PARAMETER :: GeometricFieldUserNumber=5
    96. 

  96.   INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDarcy=8
    97. 

  97.   INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumberDarcy=12
    98. 

  98.   INTEGER(CMISSIntg), PARAMETER :: ProblemUserNumber=14
    99. 

  99.   INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldUserNumberDarcy=22
    100. 

  100. 

  101.   INTEGER(CMISSIntg), PARAMETER :: ControlLoopSolidNumber=1
    102. 

  102.   INTEGER(CMISSIntg), PARAMETER :: ControlLoopFluidNumber=2
    103. 

  103.   INTEGER(CMISSIntg), PARAMETER :: ControlLoopSubiterationNumber=1
    104. 

  104.   INTEGER(CMISSIntg), PARAMETER :: SolverSolidIndex=1
    105. 

  105.   INTEGER(CMISSIntg), PARAMETER :: SolverDarcyIndex=1
    106. 

  106.   INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDarcyPorosity=1
    107. 

  107.   INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDarcyPermOverVis=2
    108. 

  108. 

  109.   INTEGER(CMISSIntg), PARAMETER :: FieldGeometryNumberOfVariables=1
    110. 

  110.   INTEGER(CMISSIntg), PARAMETER :: FieldGeometryNumberOfComponents=3
    111. 

  111. 

  112.   !Program types
    113. 

  113. 

  114.   !Program variables
    115. 

  115. 

  116.   INTEGER(CMISSIntg) :: NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS,NUMBER_GLOBAL_Z_ELEMENTS
    117. 

  117. 

  118. !   INTEGER(CMISSIntg) :: MPI_IERROR
    119. 

  119.   INTEGER(CMISSIntg) :: NumberOfComputationalNodes,NumberOfDomains,ComputationalNodeNumber
    120. 

  120. 

  121.   INTEGER(CMISSIntg) :: NUMBER_OF_DIMENSIONS
    122. 

  122. 

  123.   INTEGER(CMISSIntg) :: MAXIMUM_ITERATIONS
    124. 

  124.   INTEGER(CMISSIntg) :: RESTART_VALUE
    125. 

  125. 

  126.   INTEGER(CMISSIntg) :: EQUATIONS_DARCY_OUTPUT
    127. 

  127.   INTEGER(CMISSIntg) :: COMPONENT_NUMBER
    128. 

  128. 

  129.   INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_DARCY_OUTPUT_FREQUENCY
    130. 

  130.   INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_DARCY_OUTPUT_TYPE
    131. 

  131.   INTEGER(CMISSIntg) :: LINEAR_SOLVER_DARCY_OUTPUT_TYPE
    132. 

  132.   INTEGER(CMISSIntg) :: LINEAR_SOLVER_MAT_PROPERTIES_OUTPUT_TYPE
    133. 

  133. 

  134.   REAL(CMISSDP) :: GEOMETRY_TOLERANCE
    135. 

  135.   INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_SOLID
    136. 

  136.   REAL(CMISSDP) :: INITIAL_FIELD_DARCY(4)
    137. 

  137.   REAL(CMISSDP) :: DIVERGENCE_TOLERANCE
    138. 

  138.   REAL(CMISSDP) :: RELATIVE_TOLERANCE
    139. 

  139.   REAL(CMISSDP) :: ABSOLUTE_TOLERANCE
    140. 

  140.   REAL(CMISSDP) :: LINESEARCH_ALPHA
    141. 

  141.   REAL(CMISSDP) :: VALUE
    142. 

  142.   REAL(CMISSDP) :: POROSITY_PARAM_DARCY, PERM_OVER_VIS_PARAM_DARCY
    143. 

  143. 

  144.   LOGICAL :: EXPORT_FIELD_IO
    145. 

  145.   LOGICAL :: LINEAR_SOLVER_DARCY_DIRECT_FLAG
    146. 

  146. 

  147.   !CMISS variables
    148. 

  148. 

  149.   !Regions
    150. 

  150.   TYPE(CMISSRegionType) :: Region
    151. 

  151.   TYPE(CMISSRegionType) :: WorldRegion
    152. 

  152.   !Coordinate systems
    153. 

  153.   TYPE(CMISSCoordinateSystemType) :: CoordinateSystem
    154. 

  154.   TYPE(CMISSCoordinateSystemType) :: WorldCoordinateSystem
    155. 

  155.   !Basis
    156. 

  156. !   TYPE(CMISSBasisType) :: CubicBasis, QuadraticBasis, LinearBasis, Bases(2)
    157. 

  157.   TYPE(CMISSBasisType),allocatable :: Bases(:)
    158. 

  158.   !Meshes
    159. 

  159.   TYPE(CMISSMeshType) :: Mesh
    160. 

  160.   TYPE(CMISSGeneratedMeshType) :: GeneratedMesh
    161. 

  161. 

  162.   !Decompositions
    163. 

  163.   TYPE(CMISSDecompositionType) :: Decomposition
    164. 

  164.   !Fields
    165. 

  165.   TYPE(CMISSFieldsType) :: Fields
    166. 

  166.   !Field types
    167. 

  167.   TYPE(CMISSFieldType) :: GeometricField
    168. 

  168.   TYPE(CMISSFieldType) :: MaterialsFieldDarcy
    169. 

  169.   TYPE(CMISSFieldType) :: EquationsSetFieldDarcy
    170. 

  170.   !Boundary conditions
    171. 

  171.   TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsDarcy
    172. 

  172.   !Equations sets
    173. 

  173.   TYPE(CMISSEquationsSetType) :: EquationsSetDarcy
    174. 

  174.   !Equations
    175. 

  175.   TYPE(CMISSEquationsType) :: EquationsDarcy
    176. 

  176.   !Problems
    177. 

  177.   TYPE(CMISSProblemType) :: Problem
    178. 

  178.   !Control loops
    179. 

  179.   TYPE(CMISSControlLoopType) :: ControlLoop
    180. 

  180.   !Solvers
    181. 

  181.   TYPE(CMISSSolverType) :: DynamicSolverDarcy
    182. 

  182.   TYPE(CMISSSolverType) :: LinearSolverDarcy
    183. 

  183. !   TYPE(CMISSSolverType) :: LinearSolverSolid
    184. 

  184.   !Solver equations
    185. 

  185.   TYPE(CMISSSolverEquationsType) :: SolverEquationsDarcy
    186. 

  186.   ! nodes and elements
    187. 

  187.   TYPE(CMISSNodesType) :: Nodes
    188. 

  188.   TYPE(CMISSMeshElementsType),allocatable :: Elements(:)
    189. 

  189. 

  190.   !Other variables
    191. 

  191.   INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face1Nodes(:),Face2Nodes(:)
    192. 

  192.   INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face3Nodes(:),Face4Nodes(:)
    193. 

  193.   INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face5Nodes(:),Face6Nodes(:)
    194. 

  194.   INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face7Nodes(:),Face8Nodes(:)
    195. 

  195.   INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face9Nodes(:),Face10Nodes(:)
    196. 

  196.   INTEGER(CMISSIntg),ALLOCATABLE,TARGET :: Face11Nodes(:),Face12Nodes(:)
    197. 

  197.   INTEGER(CMISSIntg) :: FaceXi(6)
    198. 

  198.   INTEGER(CMISSIntg) :: NN,NODE,NodeDomain
    199. 

  199.   REAL(CMISSDP) :: XCoord,YCoord,ZCoord
    200. 

  200.   LOGICAL :: X_FIXED,Y_FIXED !,X_OKAY,Y_OKAY
    201. 

  201. 

  202. #ifdef WIN32
    203. 

  203.   !Quickwin type
    204. 

  204.   LOGICAL :: QUICKWIN_STATUS=.FALSE.
    205. 

  205.   TYPE(WINDOWCONFIG) :: QUICKWIN_WINDOW_CONFIG
    206. 

  206. #endif
    207. 

  207. 

  208.   !Generic CMISS variables
    209. 

  209. 

  210.   INTEGER(CMISSIntg) :: EquationsSetIndex
    211. 

  211.   INTEGER(CMISSIntg) :: Err
    212. 

  212. 

  213. 

  214.   INTEGER(CMISSIntg) :: DIAG_LEVEL_LIST(5)
    215. 

  215. !   CHARACTER(LEN=255) :: DIAG_ROUTINE_LIST(8) !,TIMING_ROUTINE_LIST(1)
    216. 

  216.   CHARACTER(LEN=255) :: DIAG_ROUTINE_LIST(1) !,TIMING_ROUTINE_LIST(1)
    217. 

  217. 

  218.   !
    219. 

  219.   !--------------------------------------------------------------------------------------------------------------------------------
    220. 

  220.   !
    221. 

  221. 

  222.   !Program variables and types (finite elasticity part)
    223. 

  223. 

  224.   !Test program parameters
    225. 

  225. 

  226.   INTEGER(CMISSIntg) :: SolidMeshComponenetNumber
    227. 

  227. 

  228.   INTEGER(CMISSIntg), PARAMETER :: FieldGeometrySolidUserNumber=1
    229. 

  229.   INTEGER(CMISSIntg), PARAMETER :: FieldGeometrySolidNumberOfVariables=1
    230. 

  230.   INTEGER(CMISSIntg), PARAMETER :: FieldGeometrySolidNumberOfComponents=3
    231. 

  231. 

  232.   INTEGER(CMISSIntg), PARAMETER :: FieldFibreSolidUserNumber=2
    233. 

  233.   INTEGER(CMISSIntg), PARAMETER :: FieldFibreSolidNumberOfVariables=1
    234. 

  234.   INTEGER(CMISSIntg), PARAMETER :: FieldFibreSolidNumberOfComponents=3
    235. 

  235. 

  236.   INTEGER(CMISSIntg), PARAMETER :: FieldMaterialSolidUserNumber=3
    237. 

  237.   INTEGER(CMISSIntg), PARAMETER :: FieldMaterialSolidNumberOfVariables=1
    238. 

  238.   INTEGER(CMISSIntg), PARAMETER :: FieldMaterialSolidNumberOfComponents=3
    239. 

  239. 

  240.   INTEGER(CMISSIntg), PARAMETER :: FieldDependentSolidUserNumber=4
    241. 

  241.   INTEGER(CMISSIntg), PARAMETER :: FieldDependentSolidNumberOfVariables=4
    242. 

  242.   INTEGER(CMISSIntg), PARAMETER :: FieldDependentSolidNumberOfComponents=4
    243. 

  243.   INTEGER(CMISSIntg), PARAMETER :: FieldDependentFluidNumberOfComponents=4  !(u,v,w,m)
    244. 

  244. 

  245.   INTEGER(CMISSIntg), PARAMETER :: EquationSetSolidUserNumber=1
    246. 

  246.   INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldSolidUserNumber=25
    247. 

  247. 

  248.   INTEGER(CMISSIntg), PARAMETER :: SolidDisplMeshComponentNumber=1
    249. 

  249.   INTEGER(CMISSIntg), PARAMETER :: SolidLagrMultMeshComponentNumber=2
    250. 

  250.   INTEGER(CMISSIntg), PARAMETER :: SolidGeometryMeshComponentNumber=SolidDisplMeshComponentNumber
    251. 

  251. 

  252.   INTEGER(CMISSIntg), PARAMETER :: DarcyVelMeshComponentNumber=SolidLagrMultMeshComponentNumber
    253. 

  253.   INTEGER(CMISSIntg), PARAMETER :: DarcyMassIncreaseMeshComponentNumber=SolidLagrMultMeshComponentNumber
    254. 

  254. !   INTEGER(CMISSIntg), PARAMETER :: DarcyGeometryMeshComponentNumber=SolidDisplMeshComponentNumber
    255. 

  255. 

  256.   INTEGER(CMISSIntg), PARAMETER :: GeneratedMeshUserNumber=32
    257. 

  257.   !Program types
    258. 

  258.   !Program variables
    259. 

  259. 

  260.   REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_START_TIME
    261. 

  261.   REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_STOP_TIME
    262. 

  262.   REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_THETA
    263. 

  263.   REAL(CMISSDP) :: DYNAMIC_SOLVER_DARCY_TIME_INCREMENT
    264. 

  264. 

  265.   !CMISS variables
    266. 

  266. 

  267.   TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsSolid
    268. 

  268.   TYPE(CMISSEquationsType) :: EquationsSolid
    269. 

  269.   TYPE(CMISSEquationsSetType) :: EquationsSetSolid
    270. 

  270.   TYPE(CMISSFieldType) :: GeometricFieldSolid,FibreFieldSolid,MaterialFieldSolid
    271. 

  271.   TYPE(CMISSFieldType) :: DependentFieldSolid,EquationsSetFieldSolid
    272. 

  272.   TYPE(CMISSSolverType) :: SolverSolid
    273. 

  273.   TYPE(CMISSSolverEquationsType) :: SolverEquationsSolid
    274. 

  274. 

  275.   !End - Program variables and types (finite elasticity part)
    276. 

  276. 

  277.   !
    278. 

  278.   !--------------------------------------------------------------------------------------------------------------------------------
    279. 

  279.   !
    280. 

  280. 

  281. 

  282. #ifdef WIN32
    283. 

  283.   !Initialise QuickWin
    284. 

  284.   QUICKWIN_WINDOW_CONFIG%TITLE="General Output" !Window title
    285. 

  285.   QUICKWIN_WINDOW_CONFIG%NUMTEXTROWS=-1 !Max possible number of rows
    286. 

  286.   QUICKWIN_WINDOW_CONFIG%MODE=QWIN$SCROLLDOWN
    287. 

  287.   !Set the window parameters
    288. 

  288.   QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
    289. 

  289.   !If attempt fails set with system estimated values
    290. 

  290.   IF(.NOT.QUICKWIN_STATUS) QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
    291. 

  291. #endif
    292. 

  292. 

  293.   !
    294. 

  294.   !================================================================================================================================
    295. 

  295.   !
    296. 

  296.   NUMBER_GLOBAL_X_ELEMENTS=1
    297. 

  297.   NUMBER_GLOBAL_Y_ELEMENTS=1
    298. 

  298.   NUMBER_GLOBAL_Z_ELEMENTS=3
    299. 

  299. 

  300.   IF(NUMBER_GLOBAL_Z_ELEMENTS==0)THEN
    301. 

  301.     NUMBER_OF_DIMENSIONS=2
    302. 

  302.   ELSE
    303. 

  303.     NUMBER_OF_DIMENSIONS=3
    304. 

  304.   ENDIF
    305. 

  305.   !PROBLEM CONTROL PANEL
    306. 

  306. 

  307. !   BASIS_XI_INTERPOLATION_SOLID=CMISS_BASIS_LINEAR_LAGRANGE_INTERPOLATION
    308. 

  308.   BASIS_XI_INTERPOLATION_SOLID=CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION
    309. 

  309.   !Set geometric tolerance
    310. 

  310.   GEOMETRY_TOLERANCE = 1.0E-12_CMISSDP
    311. 

  311.   !Set initial values
    312. 

  312.   INITIAL_FIELD_DARCY(1)=0.0_CMISSDP
    313. 

  313.   INITIAL_FIELD_DARCY(2)=0.0_CMISSDP
    314. 

  314.   INITIAL_FIELD_DARCY(3)=0.0_CMISSDP
    315. 

  315.   INITIAL_FIELD_DARCY(4)=0.0_CMISSDP
    316. 

  316.   !Set material parameters
    317. 

  317.   POROSITY_PARAM_DARCY=0.1_CMISSDP
    318. 

  318.   PERM_OVER_VIS_PARAM_DARCY=1.0_CMISSDP
    319. 

  319.   !Set output parameter
    320. 

  320.   !(NoOutput/ProgressOutput/TimingOutput/SolverOutput/SolverMatrixOutput)
    321. 

  321.   DYNAMIC_SOLVER_DARCY_OUTPUT_TYPE=CMISS_SOLVER_PROGRESS_OUTPUT
    322. 

  322.   LINEAR_SOLVER_DARCY_OUTPUT_TYPE=CMISS_SOLVER_SOLVER_OUTPUT
    323. 

  323.   !(NoOutput/TimingOutput/MatrixOutput/ElementOutput)
    324. 

  324.   EQUATIONS_DARCY_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT
    325. 

  325. 

  326.   !Set time parameter
    327. 

  327.   DYNAMIC_SOLVER_DARCY_START_TIME=0.0_CMISSDP
    328. 

  328.   DYNAMIC_SOLVER_DARCY_TIME_INCREMENT=1.0e-3_CMISSDP
    329. 

  329.   DYNAMIC_SOLVER_DARCY_STOP_TIME=2_CMISSIntg * DYNAMIC_SOLVER_DARCY_TIME_INCREMENT
    330. 

  330.   DYNAMIC_SOLVER_DARCY_THETA=1.0_CMISSDP !2.0_CMISSDP/3.0_CMISSDP
    331. 

  331.   !Set result output parameter
    332. 

  332.   DYNAMIC_SOLVER_DARCY_OUTPUT_FREQUENCY=1
    333. 

  333.   !Set solver parameters
    334. 

  334.   LINEAR_SOLVER_DARCY_DIRECT_FLAG=.TRUE.
    335. 

  335.   RELATIVE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-05_CMISSDP
    336. 

  336.   ABSOLUTE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-10_CMISSDP
    337. 

  337.   DIVERGENCE_TOLERANCE=1.0E5_CMISSDP !default: 1.0E5
    338. 

  338.   MAXIMUM_ITERATIONS=10000_CMISSIntg !default: 100000
    339. 

  339.   RESTART_VALUE=30_CMISSIntg !default: 30
    340. 

  340.   LINESEARCH_ALPHA=1.0_CMISSDP
    341. 

  341. 

  342. 

  343.   !
    344. 

  344.   !================================================================================================================================
    345. 

  345.   !
    346. 

  346. 

  347.   !INITIALISE OPENCMISS
    348. 

  348. 

  349.   CALL CMISSInitialise(WorldCoordinateSystem,WorldRegion,Err)
    350. 

  350. 

  351.   CALL CMISSErrorHandlingModeSet(CMISS_ERRORS_TRAP_ERROR,Err)
    352. 

  352. 

  353.   !
    354. 

  354.   !================================================================================================================================
    355. 

  355.   !
    356. 

  356. 

  357.   !Set diagnostics
    358. 

  358. 

  359.   DIAG_LEVEL_LIST(1)=1
    360. 

  360.   DIAG_LEVEL_LIST(2)=2
    361. 

  361.   DIAG_LEVEL_LIST(3)=3
    362. 

  362.   DIAG_LEVEL_LIST(4)=4
    363. 

  363.   DIAG_LEVEL_LIST(5)=5
    364. 

  364. 

  365.   DIAG_ROUTINE_LIST(1)="WRITE_IP_INFO"
    366. 

  366. !   DIAG_ROUTINE_LIST(2)="FINITE_ELASTICITY_GAUSS_CAUCHY_TENSOR"
    367. 

  367. 

  368.   !CMISS_ALL_DIAG_TYPE/CMISS_IN_DIAG_TYPE/CMISS_FROM_DIAG_TYPE
    369. 

  369.   CALL CMISSDiagnosticsSetOn(CMISS_IN_DIAG_TYPE,DIAG_LEVEL_LIST,"Diagnostics",DIAG_ROUTINE_LIST,Err)
    370. 

  370. 

  371.   !CMISS_ALL_TIMING_TYPE/CMISS_IN_TIMING_TYPE/CMISS_FROM_TIMING_TYPE
    372. 

  372.   !TIMING_ROUTINE_LIST(1)="PROBLEM_FINITE_ELEMENT_CALCULATE"
    373. 

  373.   !CALL TIMING_SET_ON(IN_TIMING_TYPE,.TRUE.,"",TIMING_ROUTINE_LIST,ERR,ERROR,*999)
    374. 

  374. 

  375.   !
    376. 

  376.   !================================================================================================================================
    377. 

  377.   !
    378. 

  378. 

  379.   !Get the number of computational nodes and this computational node number
    380. 

  380.   CALL CMISSComputationalNumberOfNodesGet(NumberOfComputationalNodes,Err)
    381. 

  381.   CALL CMISSComputationalNodeNumberGet(ComputationalNodeNumber,Err)
    382. 

  382. 

  383.   NumberOfDomains = NumberOfComputationalNodes
    384. 

  384.   write(*,*) "NumberOfDomains = ",NumberOfDomains
    385. 

  385. 

  386.   !
    387. 

  387.   !================================================================================================================================
    388. 

  388.   !
    389. 

  389. 

  390.   !COORDINATE SYSTEM
    391. 

  391. 

  392.   CALL CMISSCoordinateSystem_Initialise(CoordinateSystem,Err)
    393. 

  393.   CALL CMISSCoordinateSystem_CreateStart(CoordinateSystemUserNumber,CoordinateSystem,Err)
    394. 

  394.   CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,NUMBER_OF_DIMENSIONS,Err)
    395. 

  395.   CALL CMISSCoordinateSystem_CreateFinish(CoordinateSystem,Err)
    396. 

  396. 

  397.   !
    398. 

  398.   !================================================================================================================================
    399. 

  399.   !
    400. 

  400. 

  401.   !REGION
    402. 

  402.   !For a volume-coupled problem, solid and fluid are based in the same region
    403. 

  403. 

  404.   CALL CMISSRegion_Initialise(Region,Err)
    405. 

  405.   CALL CMISSRegion_CreateStart(RegionUserNumber,WorldRegion,Region,Err)
    406. 

  406.   CALL CMISSRegion_CoordinateSystemSet(Region,CoordinateSystem,Err)
    407. 

  407.   CALL CMISSRegion_CreateFinish(Region,Err)
    408. 

  408. 

  409.   !
    410. 

  410.   !================================================================================================================================
    411. 

  411.   !
    412. 

  412. 

  413. call READ_MESH('input/example-mesh-3els',MeshUserNumber,Region, Mesh,Bases,Nodes,Elements)
    414. 

  414. !   !BASES
    415. 

  415. !   !Define basis functions
    416. 

  416. !   CALL CMISSBasis_Initialise(LinearBasis,Err)
    417. 

  417. !   CALL CMISSBasis_CreateStart(LinearBasisUserNumber,LinearBasis,Err)
    418. 

  418. !   CALL CMISSBasis_QuadratureNumberOfGaussXiSet(LinearBasis, &
    419. 

  419. !     & (/CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME/),Err)
    420. 

  420. !   !CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(LinearBasis,.true.,Err)
    421. 

  421. !   CALL CMISSBasis_CreateFinish(LinearBasis,Err)
    422. 

  422. ! 
    423. 

  423. !   CALL CMISSBasis_Initialise(QuadraticBasis,Err)
    424. 

  424. !   CALL CMISSBasis_CreateStart(QuadraticBasisUserNumber,QuadraticBasis,Err)
    425. 

  425. !   CALL CMISSBasis_InterpolationXiSet(QuadraticBasis,(/CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION, &
    426. 

  426. !     & CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION,CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION/),Err)
    427. 

  427. !   CALL CMISSBasis_QuadratureNumberOfGaussXiSet(QuadraticBasis, &
    428. 

  428. !     & (/CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME/),Err)
    429. 

  429. !   !CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(QuadraticBasis,.true.,Err)
    430. 

  430. !   CALL CMISSBasis_CreateFinish(QuadraticBasis,Err)
    431. 

  431. ! 
    432. 

  432. !   CALL CMISSBasis_Initialise(CubicBasis,Err)
    433. 

  433. !   CALL CMISSBasis_CreateStart(CubicBasisUserNumber,CubicBasis,Err)
    434. 

  434. !   CALL CMISSBasis_InterpolationXiSet(CubicBasis,(/CMISS_BASIS_CUBIC_LAGRANGE_INTERPOLATION, &
    435. 

  435. !     & CMISS_BASIS_CUBIC_LAGRANGE_INTERPOLATION,CMISS_BASIS_CUBIC_LAGRANGE_INTERPOLATION/),Err)
    436. 

  436. !   CALL CMISSBasis_QuadratureNumberOfGaussXiSet(CubicBasis, &
    437. 

  437. !     & (/CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME,CMISS_BASIS_HIGH_QUADRATURE_SCHEME/),Err)
    438. 

  438. !   !CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(CubicBasis,.true.,Err) !Enable 3D interpolation on faces
    439. 

  439. !   CALL CMISSBasis_CreateFinish(CubicBasis,Err)
    440. 

  440. ! 
    441. 

  441. !   !LinearBasis/QuadraticBasis/CubicBasis
    442. 

  442. !   Bases(1)=QuadraticBasis
    443. 

  443. !   Bases(2)=LinearBasis
    444. 

  444. ! 
    445. 

  445. !   !Start the creation of a generated mesh in the region
    446. 

  446. !   CALL CMISSGeneratedMesh_Initialise(GeneratedMesh,Err)
    447. 

  447. !   CALL CMISSGeneratedMesh_CreateStart(GeneratedMeshUserNumber,Region,GeneratedMesh,Err)
    448. 

  448. !   CALL CMISSGeneratedMesh_TypeSet(GeneratedMesh,CMISS_GENERATED_MESH_REGULAR_MESH_TYPE,Err)
    449. 

  449. !   CALL CMISSGeneratedMesh_BasisSet(GeneratedMesh,Bases,Err)
    450. 

  450. !   IF(NUMBER_GLOBAL_Z_ELEMENTS==0) THEN
    451. 

  451. !     CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/X_DIM,Y_DIM/),Err)
    452. 

  452. !     CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS/),Err)
    453. 

  453. !   ELSE
    454. 

  454. !     CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/X_DIM,Y_DIM,Z_DIM/),Err)
    455. 

  455. !     CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS, &
    456. 

  456. !       & NUMBER_GLOBAL_Z_ELEMENTS/),Err)
    457. 

  457. !   ENDIF
    458. 

  458. !   CALL CMISSMesh_Initialise(Mesh,Err)
    459. 

  459. !   CALL CMISSGeneratedMesh_CreateFinish(GeneratedMesh,MeshUserNumber,Mesh,Err)
    460. 

  460. ! -------------------------------------------------------------------------------
    461. 

  461.   
    462. 

  462. 

  463.   !GEOMETRIC FIELD
    464. 

  464. 

  465.   !Create a decomposition:
    466. 

  466.   CALL CMISSDecomposition_Initialise(Decomposition,Err)
    467. 

  467.   CALL CMISSDecomposition_CreateStart(DecompositionUserNumber,Mesh,Decomposition,Err)
    468. 

  468.   !Set the decomposition to be a general decomposition with the specified number of domains
    469. 

  469.   CALL CMISSDecomposition_TypeSet(Decomposition,CMISS_DECOMPOSITION_CALCULATED_TYPE,Err)
    470. 

  470.   CALL CMISSDecomposition_NumberOfDomainsSet(Decomposition,NumberOfDomains,Err)
    471. 

  471.   CALL CMISSDecomposition_CreateFinish(Decomposition,Err)
    472. 

  472. 

  473.   CALL CMISSField_Initialise(GeometricField,Err)
    474. 

  474.   CALL CMISSField_CreateStart(GeometricFieldUserNumber,Region,GeometricField,Err)
    475. 

  475.   CALL CMISSField_MeshDecompositionSet(GeometricField,Decomposition,Err)
    476. 

  476.   CALL CMISSField_TypeSet(GeometricField,CMISS_FIELD_GEOMETRIC_TYPE,Err)
    477. 

  477.   CALL CMISSField_NumberOfVariablesSet(GeometricField,FieldGeometryNumberOfVariables,Err)
    478. 

  478.   CALL CMISSField_NumberOfComponentsSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,3,Err)  
    479. 

  479.   CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidGeometryMeshComponentNumber,Err)
    480. 

  480.   CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidGeometryMeshComponentNumber,Err)
    481. 

  481.   CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidGeometryMeshComponentNumber,Err)
    482. 

  482.   CALL CMISSField_CreateFinish(GeometricField,Err)
    483. 

  483. CALL READ_NODES('input/example-nodes-3els',GeometricField)
    484. 

  484. !   CALL CMISSGeneratedMesh_GeometricParametersCalculate(GeneratedMesh,GeometricField,Err)
    485. 

  485. 

  486.   !--------------------------------------------------------------------------------------------------------------------------------
    487. 

  487.   ! Solid
    488. 

  488. 

  489.   !Create a decomposition
    490. 

  490. 

  491.   !Create a field to put the geometry (defualt is geometry)
    492. 

  492. 

  493.   SolidMeshComponenetNumber = SolidGeometryMeshComponentNumber
    494. 

  494. 

  495.   CALL CMISSField_Initialise(GeometricFieldSolid,Err)
    496. 

  496.   CALL CMISSField_CreateStart(FieldGeometrySolidUserNumber,Region,GeometricFieldSolid,Err)
    497. 

  497.   CALL CMISSField_MeshDecompositionSet(GeometricFieldSolid,Decomposition,Err)
    498. 

  498.   CALL CMISSField_TypeSet(GeometricFieldSolid,CMISS_FIELD_GEOMETRIC_TYPE,Err)
    499. 

  499.   CALL CMISSField_NumberOfVariablesSet(GeometricFieldSolid,FieldGeometrySolidNumberOfVariables,Err)
    500. 

  500.   CALL CMISSField_NumberOfComponentsSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldGeometrySolidNumberOfComponents,Err)
    501. 

  501.   CALL CMISSField_ComponentMeshComponentSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidMeshComponenetNumber,Err)
    502. 

  502.   CALL CMISSField_ComponentMeshComponentSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidMeshComponenetNumber,Err)
    503. 

  503.   CALL CMISSField_ComponentMeshComponentSet(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidMeshComponenetNumber,Err)
    504. 

  504.   CALL CMISSField_CreateFinish(GeometricFieldSolid,Err)
    505. 

  505.   !Set the mesh component to be used by the field components.
    506. 

  506. CALL READ_NODES('input/example-nodes-3els',GeometricFieldSolid)
    507. 

  507. !  CALL CMISSGeneratedMesh_GeometricParametersCalculate(GeneratedMesh,GeometricFieldSolid,Err)
    508. 

  508. 

  509.   !Create a fibre field and attach it to the geometric field
    510. 

  510.   CALL CMISSField_Initialise(FibreFieldSolid,Err)
    511. 

  511.   CALL CMISSField_CreateStart(FieldFibreSolidUserNumber,Region,FibreFieldSolid,Err)
    512. 

  512.   CALL CMISSField_TypeSet(FibreFieldSolid,CMISS_FIELD_FIBRE_TYPE,Err)
    513. 

  513.   CALL CMISSField_MeshDecompositionSet(FibreFieldSolid,Decomposition,Err)
    514. 

  514.   CALL CMISSField_GeometricFieldSet(FibreFieldSolid,GeometricFieldSolid,Err)
    515. 

  515.   CALL CMISSField_NumberOfVariablesSet(FibreFieldSolid,FieldFibreSolidNumberOfVariables,Err)
    516. 

  516.   CALL CMISSField_NumberOfComponentsSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldFibreSolidNumberOfComponents,Err)
    517. 

  517.   CALL CMISSField_ComponentMeshComponentSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidGeometryMeshComponentNumber,Err)
    518. 

  518.   CALL CMISSField_ComponentMeshComponentSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidGeometryMeshComponentNumber,Err)
    519. 

  519.   CALL CMISSField_ComponentMeshComponentSet(FibreFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidGeometryMeshComponentNumber,Err)
    520. 

  520.   CALL CMISSField_CreateFinish(FibreFieldSolid,Err)
    521. 

  521. 

  522.   ! end Solid
    523. 

  523.   !--------------------------------------------------------------------------------------------------------------------------------
    524. 

  524. 

  525.   !
    526. 

  526.   !================================================================================================================================
    527. 

  527.   !
    528. 

  528. 

  529.   !EQUATIONS SETS
    530. 

  530. 

  531.   !Create the equations set for ALE Darcy
    532. 

  532.   CALL CMISSField_Initialise(EquationsSetFieldDarcy,Err)
    533. 

  533.   CALL CMISSEquationsSet_Initialise(EquationsSetDarcy,Err)
    534. 

  534.   CALL CMISSEquationsSet_CreateStart(EquationsSetUserNumberDarcy,Region,GeometricField,CMISS_EQUATIONS_SET_FLUID_MECHANICS_CLASS, &
    535. 

  535.     & CMISS_EQUATIONS_SET_DARCY_EQUATION_TYPE,CMISS_EQUATIONS_SET_INCOMPRESS_ELASTICITY_DRIVEN_DARCY_SUBTYPE,&
    536. 

  536.     & EquationsSetFieldUserNumberDarcy,EquationsSetFieldDarcy,EquationsSetDarcy,Err)
    537. 

  537.   CALL CMISSEquationsSet_CreateFinish(EquationsSetDarcy,Err)
    538. 

  538. 

  539.   !Create the equations set for the solid
    540. 

  540.   CALL CMISSField_Initialise(EquationsSetFieldSolid,Err)
    541. 

  541.   CALL CMISSEquationsSet_Initialise(EquationsSetSolid,Err)
    542. 

  542.   CALL CMISSEquationsSet_CreateStart(EquationSetSolidUserNumber,Region,FibreFieldSolid,CMISS_EQUATIONS_SET_ELASTICITY_CLASS, &
    543. 

  543.     & CMISS_EQUATIONS_SET_FINITE_ELASTICITY_TYPE,CMISS_EQUATIONS_SET_INCOMPRESS_ELASTICITY_DRIVEN_DARCY_SUBTYPE,&
    544. 

  544.     & EquationsSetFieldSolidUserNumber,EquationsSetFieldSolid,EquationsSetSolid,Err)
    545. 

  545.   CALL CMISSEquationsSet_CreateFinish(EquationsSetSolid,Err)
    546. 

  546. 

  547.   !--------------------------------------------------------------------------------------------------------------------------------
    548. 

  548.   ! Solid Materials Field
    549. 

  549. 

  550.   !Create a material field and attach it to the geometric field
    551. 

  551.   CALL CMISSField_Initialise(MaterialFieldSolid,Err)
    552. 

  552.   !
    553. 

  553.   CALL CMISSField_CreateStart(FieldMaterialSolidUserNumber,Region,MaterialFieldSolid,Err)
    554. 

  554.   !
    555. 

  555.   CALL CMISSField_TypeSet(MaterialFieldSolid,CMISS_FIELD_MATERIAL_TYPE,Err)
    556. 

  556.   CALL CMISSField_MeshDecompositionSet(MaterialFieldSolid,Decomposition,Err)
    557. 

  557.   CALL CMISSField_GeometricFieldSet(MaterialFieldSolid,GeometricFieldSolid,Err)
    558. 

  558.   CALL CMISSField_NumberOfVariablesSet(MaterialFieldSolid,FieldMaterialSolidNumberOfVariables,Err)
    559. 

  559.   CALL CMISSField_NumberOfComponentsSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldMaterialSolidNumberOfComponents,Err)
    560. 

  560.   CALL CMISSField_ComponentMeshComponentSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidGeometryMeshComponentNumber,Err)
    561. 

  561.   CALL CMISSField_ComponentMeshComponentSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidGeometryMeshComponentNumber,Err)
    562. 

  562.   CALL CMISSField_ComponentMeshComponentSet(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidGeometryMeshComponentNumber,Err)
    563. 

  563.   !
    564. 

  564.   CALL CMISSField_CreateFinish(MaterialFieldSolid,Err)
    565. 

  565. 

  566.   !Set material parameters
    567. 

  567.   CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1, &
    568. 

  568.     & 2.0_CMISSDP,Err)
    569. 

  569. !   CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,2.0e3_CMISSDP,Err)
    570. 

  570.   CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,2, &
    571. 

  571.     & 6.0_CMISSDP,Err)
    572. 

  572. !   CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,2,33.0_CMISSDP,Err)
    573. 

  573.   CALL CMISSField_ComponentValuesInitialise(MaterialFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,3, &
    574. 

  574.     & 10.0_CMISSDP,Err)
    575. 

  575. 

  576. 

  577.   CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetSolid,FieldMaterialSolidUserNumber,MaterialFieldSolid,Err)
    578. 

  578.   CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetSolid,Err)
    579. 

  579. 

  580.   ! end Solid
    581. 

  581.   !--------------------------------------------------------------------------------------------------------------------------------
    582. 

  582. 

  583. 

  584.   !
    585. 

  585.   !================================================================================================================================
    586. 

  586.   !
    587. 

  587. 

  588.   !DEPENDENT FIELDS
    589. 

  589. 

  590.   !--------------------------------------------------------------------------------------------------------------------------------
    591. 

  591.   ! Solid
    592. 

  592. 

  593.   !Create a dependent field with four variables (U, DelUDelN = solid, V, DelVDelN = Darcy) and four components
    594. 

  594.   CALL CMISSField_Initialise(DependentFieldSolid,Err)
    595. 

  595.   !
    596. 

  596.   CALL CMISSField_CreateStart(FieldDependentSolidUserNumber,Region,DependentFieldSolid,Err)
    597. 

  597.   !
    598. 

  598.   CALL CMISSField_TypeSet(DependentFieldSolid,CMISS_FIELD_GENERAL_TYPE,Err)
    599. 

  599.   CALL CMISSField_MeshDecompositionSet(DependentFieldSolid,Decomposition,Err)
    600. 

  600.   CALL CMISSField_GeometricFieldSet(DependentFieldSolid,GeometricFieldSolid,Err)
    601. 

  601.   CALL CMISSField_DependentTypeSet(DependentFieldSolid,CMISS_FIELD_DEPENDENT_TYPE,Err)
    602. 

  602.   CALL CMISSField_NumberOfVariablesSet(DependentFieldSolid,FieldDependentSolidNumberOfVariables,Err)
    603. 

  603.   CALL CMISSField_VariableTypesSet(DependentFieldSolid,(/CMISS_FIELD_U_VARIABLE_TYPE, &
    604. 

  604.     & CMISS_FIELD_DELUDELN_VARIABLE_TYPE,CMISS_FIELD_V_VARIABLE_TYPE,CMISS_FIELD_DELVDELN_VARIABLE_TYPE/),Err)
    605. 

  605.   CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,FieldDependentSolidNumberOfComponents,Err)
    606. 

  606.   CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE, &
    607. 

  607.     & FieldDependentSolidNumberOfComponents,Err)
    608. 

  608.   CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,FieldDependentFluidNumberOfComponents,Err)
    609. 

  609.   CALL CMISSField_NumberOfComponentsSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE, &
    610. 

  610.     & FieldDependentFluidNumberOfComponents,Err)
    611. 

  611.   !
    612. 

  612.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,SolidDisplMeshComponentNumber,Err)
    613. 

  613.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,2,SolidDisplMeshComponentNumber,Err)
    614. 

  614.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,3,SolidDisplMeshComponentNumber,Err)
    615. 

  615.   CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4, &
    616. 

  616.     & CMISS_FIELD_NODE_BASED_INTERPOLATION, &
    617. 

  617.     & Err)
    618. 

  618. !   CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4,CMISS_FIELD_ELEMENT_BASED_INTERPOLATION,Err)
    619. 

  619. !   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4,SolidMeshComponenetNumber,Err)
    620. 

  620.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,4,SolidLagrMultMeshComponentNumber,Err)
    621. 

  621.   !
    622. 

  622.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,1, &
    623. 

  623.     & SolidDisplMeshComponentNumber, &
    624. 

  624.     & Err)
    625. 

  625.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,2, &
    626. 

  626.     & SolidDisplMeshComponentNumber, &
    627. 

  627.     & Err)
    628. 

  628.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,3, &
    629. 

  629.     & SolidDisplMeshComponentNumber, &
    630. 

  630.     & Err)
    631. 

  631.   CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4, &
    632. 

  632.     & CMISS_FIELD_NODE_BASED_INTERPOLATION,Err)
    633. 

  633. !   CALL CMISSField_ComponentInterpolationSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4, &
    634. 

  634. !     & CMISS_FIELD_ELEMENT_BASED_INTERPOLATION,Err)
    635. 

  635. !   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4,SolidMeshComponenetNumber,Err)
    636. 

  636.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,4, &
    637. 

  637.     & SolidLagrMultMeshComponentNumber, &
    638. 

  638.     & Err)
    639. 

  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. 

  641.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,DarcyVelMeshComponentNumber,Err)
    642. 

  642.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,2,DarcyVelMeshComponentNumber,Err)
    643. 

  643.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,3,DarcyVelMeshComponentNumber,Err)
    644. 

  644. !   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,4,MESH_COMPONENT_NUMBER_PRESSURE,Err)
    645. 

  645.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,4, &
    646. 

  646.     & DarcyMassIncreaseMeshComponentNumber, &
    647. 

  647.     & Err)
    648. 

  648.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,1,DarcyVelMeshComponentNumber, &
    649. 

  649.     & Err)
    650. 

  650.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,2,DarcyVelMeshComponentNumber, &
    651. 

  651.     & Err)
    652. 

  652.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,3,DarcyVelMeshComponentNumber, &
    653. 

  653.     & Err)
    654. 

  654. !   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,4,MESH_COMPONENT_NUMBER_PRESSURE,Err)
    655. 

  655.   CALL CMISSField_ComponentMeshComponentSet(DependentFieldSolid,CMISS_FIELD_DELVDELN_VARIABLE_TYPE,4, &
    656. 

  656.     & DarcyMassIncreaseMeshComponentNumber,Err)
    657. 

  657. 

  658.   CALL CMISSField_CreateFinish(DependentFieldSolid,Err)
    659. 

  659.   !
    660. 

  660.   CALL CMISSEquationsSet_DependentCreateStart(EquationsSetSolid,FieldDependentSolidUserNumber,DependentFieldSolid,Err)
    661. 

  661.   CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetSolid,Err)
    662. 

  662. 

  663.   ! end Solid
    664. 

  664.   !--------------------------------------------------------------------------------------------------------------------------------
    665. 

  665. 

  666. 

  667.   !Create the equations set dependent field variables for ALE Darcy
    668. 

  668.   CALL CMISSEquationsSet_DependentCreateStart(EquationsSetDarcy,FieldDependentSolidUserNumber, & ! ??? UserNumber ???
    669. 

  669.     & DependentFieldSolid,Err)
    670. 

  670.   CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetDarcy,Err)
    671. 

  671. 

  672.   !Initialise dependent field (velocity components,mass increase)
    673. 

  673.   DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS+1
    674. 

  674.     CALL CMISSField_ComponentValuesInitialise(DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
    675. 

  675.       & COMPONENT_NUMBER,INITIAL_FIELD_DARCY(COMPONENT_NUMBER),Err)
    676. 

  676.   ENDDO
    677. 

  677. 

  678. 

  679.   !
    680. 

  680.   !================================================================================================================================
    681. 

  681.   !
    682. 

  682. 

  683. call READ_FIELD('input/example-field-darcy',MaterialsFieldUserNumberDarcy,Region,GeometricField,MaterialsFieldDarcy)
    684. 

  684. CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetDarcy,MaterialsFieldUserNumberDarcy,MaterialsFieldDarcy,Err)
    685. 

  685. CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetDarcy,Err)
    686. 

  686. !   !MATERIALS FIELDS
    687. 

  687. ! 
    688. 

  688. !   !Create the equations set materials field variables for ALE Darcy
    689. 

  689. !   CALL CMISSField_Initialise(MaterialsFieldDarcy,Err)
    690. 

  690. !   CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetDarcy,MaterialsFieldUserNumberDarcy, &
    691. 

  691. !     & MaterialsFieldDarcy,Err)
    692. 

  692. !   !Finish the equations set materials field variables
    693. 

  693. !   CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetDarcy,Err)
    694. 

  694. !   CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDarcy,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
    695. 

  695. !     & MaterialsFieldUserNumberDarcyPorosity,POROSITY_PARAM_DARCY,Err)
    696. 

  696. !   CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDarcy,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
    697. 

  697. !     & MaterialsFieldUserNumberDarcyPermOverVis,PERM_OVER_VIS_PARAM_DARCY,Err)
    698. 

  698. 

  699.   !
    700. 

  700.   !================================================================================================================================
    701. 

  701.   !
    702. 

  702. 

  703.   !EQUATIONS SET EQUATIONS
    704. 

  704. 

  705.   !Darcy
    706. 

  706.   CALL CMISSEquations_Initialise(EquationsDarcy,Err)
    707. 

  707.   CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetDarcy,EquationsDarcy,Err)
    708. 

  708.   CALL CMISSEquations_SparsityTypeSet(EquationsDarcy,CMISS_EQUATIONS_SPARSE_MATRICES,Err)
    709. 

  709.   CALL CMISSEquations_OutputTypeSet(EquationsDarcy,EQUATIONS_DARCY_OUTPUT,Err)
    710. 

  710.   CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetDarcy,Err)
    711. 

  711. 

  712.   !Solid
    713. 

  713.   CALL CMISSEquations_Initialise(EquationsSolid,Err)
    714. 

  714.   CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetSolid,EquationsSolid,Err)
    715. 

  715.   CALL CMISSEquations_SparsityTypeSet(EquationsSolid,CMISS_EQUATIONS_SPARSE_MATRICES,Err)
    716. 

  716.   CALL CMISSEquations_OutputTypeSet(EquationsSolid,CMISS_EQUATIONS_NO_OUTPUT,Err)
    717. 

  717.   CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetSolid,Err)
    718. 

  718. 

  719.   !
    720. 

  720.   !================================================================================================================================
    721. 

  721.   !
    722. 

  722. 

  723.   !--------------------------------------------------------------------------------------------------------------------------------
    724. 

  724.   ! Solid
    725. 

  725. 

  726.   !Initialise dependent field from undeformed geometry and displacement bcs and set hydrostatic pressure
    727. 

  727.   CALL CMISSField_ParametersToFieldParametersComponentCopy(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE, &
    728. 

  728.     & CMISS_FIELD_VALUES_SET_TYPE, &
    729. 

  729.     & 1,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,Err)
    730. 

  730.   CALL CMISSField_ParametersToFieldParametersComponentCopy(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE, &
    731. 

  731.     & CMISS_FIELD_VALUES_SET_TYPE, &
    732. 

  732.     & 2,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,2,Err)
    733. 

  733.   CALL CMISSField_ParametersToFieldParametersComponentCopy(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE, &
    734. 

  734.     & CMISS_FIELD_VALUES_SET_TYPE, &
    735. 

  735.     & 3,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,3,Err)
    736. 

  736.   CALL CMISSField_ComponentValuesInitialise(DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,4, &
    737. 

  737.     & 0.0_CMISSDP, &
    738. 

  738.     & Err)
    739. 

  739. 

  740. 

  741.   !
    742. 

  742.   !================================================================================================================================
    743. 

  743.   !
    744. 

  744. 

  745.   !PROBLEMS
    746. 

  746. 

  747.   CALL CMISSProblem_Initialise(Problem,Err)
    748. 

  748.   CALL CMISSControlLoop_Initialise(ControlLoop,Err)
    749. 

  749.   CALL CMISSProblem_CreateStart(ProblemUserNumber,Problem,Err)
    750. 

  750.   CALL CMISSProblem_SpecificationSet(Problem,CMISS_PROBLEM_MULTI_PHYSICS_CLASS,CMISS_PROBLEM_FINITE_ELASTICITY_DARCY_TYPE, &
    751. 

  751.     & CMISS_PROBLEM_QUASISTATIC_ELASTICITY_TRANSIENT_DARCY_SUBTYPE,Err)
    752. 

  752.   CALL CMISSProblem_CreateFinish(Problem,Err)
    753. 

  753. 

  754.   CALL CMISSProblem_ControlLoopCreateStart(Problem,Err)
    755. 

  755.   CALL CMISSProblem_ControlLoopGet(Problem,CMISS_CONTROL_LOOP_NODE,ControlLoop,Err)
    756. 

  756. !   CALL CMISSControlLoop_MaximumIterationsSet(ControlLoop,2,Err)
    757. 

  757.   CALL CMISSControlLoop_TimesSet(ControlLoop,DYNAMIC_SOLVER_DARCY_START_TIME,DYNAMIC_SOLVER_DARCY_STOP_TIME, &
    758. 

  758.     & DYNAMIC_SOLVER_DARCY_TIME_INCREMENT,Err)
    759. 

  759.   CALL CMISSControlLoop_TimeOutputSet(ControlLoop,DYNAMIC_SOLVER_DARCY_OUTPUT_FREQUENCY,Err)
    760. 

  760. !   CALL CMISSControlLoop_OutputTypeSet(ControlLoop,CMISS_CONTROL_LOOP_PROGRESS_OUTPUT,Err)
    761. 

  761.   CALL CMISSProblem_ControlLoopCreateFinish(Problem,Err)
    762. 

  762. 

  763.   !
    764. 

  764.   !================================================================================================================================
    765. 

  765.   !
    766. 

  766. 

  767.   !SOLVERS
    768. 

  768. 

  769.   CALL CMISSSolver_Initialise(SolverSolid,Err)
    770. 

  770.   CALL CMISSSolver_Initialise(DynamicSolverDarcy,Err)
    771. 

  771.   CALL CMISSSolver_Initialise(LinearSolverDarcy,Err)
    772. 

  772. 

  773.   CALL CMISSProblem_SolversCreateStart(Problem,Err)
    774. 

  774. 

  775.   ! Solid
    776. 

  776.   CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopSolidNumber,CMISS_CONTROL_LOOP_NODE/), &
    777. 

  777.     & SolverSolidIndex,SolverSolid,Err)
    778. 

  778.   CALL CMISSSolver_OutputTypeSet(SolverSolid,CMISS_SOLVER_PROGRESS_OUTPUT,Err)
    779. 

  779. !   CALL CMISSSolver_NewtonJacobianCalculationTypeSet(SolverSolid,CMISS_SOLVER_NEWTON_JACOBIAN_FD_CALCULATED,Err)
    780. 

  780.   CALL CMISSSolver_NewtonJacobianCalculationTypeSet(SolverSolid,CMISS_SOLVER_NEWTON_JACOBIAN_EQUATIONS_CALCULATED,Err)
    781. 

  781. 

  782.   CALL CMISSSolver_NewtonAbsoluteToleranceSet(SolverSolid,ABSOLUTE_TOLERANCE,Err)
    783. 

  783.   CALL CMISSSolver_NewtonRelativeToleranceSet(SolverSolid,RELATIVE_TOLERANCE,Err)
    784. 

  784.   CALL CMISSSolver_NewtonMaximumIterationsSet(SolverSolid,MAXIMUM_ITERATIONS,Err)
    785. 

  785. 

  786. !   CALL CMISSSolverNonLinearTypeSet(SolverSolid,CMISS_SOLVER_NONLINEAR_NEWTON,Err)
    787. 

  787. !   CALL CMISSSolver_LibraryTypeSet(SolverSolid,CMISS_SOLVER_PETSC_LIBRARY,Err)
    788. 

  788. 

  789. !   CALL CMISSSolver_NewtonLinearSolverGet(SolverSolid,LinearSolverSolid,Err)
    790. 

  790. !   CALL CMISSSolver_LinearTypeSet(LinearSolverSolid,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err)
    791. 

  791. 

  792. 

  793.   !Darcy
    794. 

  794.   CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopFluidNumber,CMISS_CONTROL_LOOP_NODE/), &
    795. 

  795.     & SolverDarcyIndex,DynamicSolverDarcy,Err)
    796. 

  796.   CALL CMISSSolver_OutputTypeSet(DynamicSolverDarcy,DYNAMIC_SOLVER_DARCY_OUTPUT_TYPE,Err)
    797. 

  797.   CALL CMISSSolver_DynamicThetaSet(DynamicSolverDarcy,DYNAMIC_SOLVER_DARCY_THETA,Err)
    798. 

  798. !   CALL CMISSSolverDynamicDynamicSet(DynamicSolverDarcy,.TRUE.,Err)
    799. 

  799.   CALL CMISSSolver_DynamicLinearSolverGet(DynamicSolverDarcy,LinearSolverDarcy,Err)
    800. 

  800.   IF(LINEAR_SOLVER_DARCY_DIRECT_FLAG) THEN
    801. 

  801.     CALL CMISSSolver_LinearTypeSet(LinearSolverDarcy,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err)
    802. 

  802.     CALL CMISSSolver_LibraryTypeSet(LinearSolverDarcy,CMISS_SOLVER_MUMPS_LIBRARY,Err)
    803. 

  803.   ELSE
    804. 

  804.     CALL CMISSSolver_LinearTypeSet(LinearSolverDarcy,CMISS_SOLVER_LINEAR_ITERATIVE_SOLVE_TYPE,Err)
    805. 

  805.     CALL CMISSSolver_LinearIterativeMaximumIterationsSet(LinearSolverDarcy,MAXIMUM_ITERATIONS,Err)
    806. 

  806.     CALL CMISSSolver_LinearIterativeDivergenceToleranceSet(LinearSolverDarcy,DIVERGENCE_TOLERANCE,Err)
    807. 

  807.     CALL CMISSSolver_LinearIterativeRelativeToleranceSet(LinearSolverDarcy,RELATIVE_TOLERANCE,Err)
    808. 

  808.     CALL CMISSSolver_LinearIterativeAbsoluteToleranceSet(LinearSolverDarcy,ABSOLUTE_TOLERANCE,Err)
    809. 

  809.     CALL CMISSSolver_LinearIterativeGMRESRestartSet(LinearSolverDarcy,RESTART_VALUE,Err)
    810. 

  810.   ENDIF
    811. 

  811. 

  812.   CALL CMISSProblem_SolversCreateFinish(Problem,Err)
    813. 

  813. 

  814.   !
    815. 

  815.   !================================================================================================================================
    816. 

  816.   !
    817. 

  817. 

  818.   !SOLVER EQUATIONS
    819. 

  819. 

  820.   CALL CMISSSolver_Initialise(SolverSolid,Err)
    821. 

  821.   CALL CMISSSolver_Initialise(LinearSolverDarcy,Err)
    822. 

  822. 

  823.   CALL CMISSSolverEquations_Initialise(SolverEquationsSolid,Err)
    824. 

  824.   CALL CMISSSolverEquations_Initialise(SolverEquationsDarcy,Err)
    825. 

  825. 

  826.   CALL CMISSProblem_SolverEquationsCreateStart(Problem,Err)
    827. 

  827.   !
    828. 

  828.   !Get the finite elasticity solver equations
    829. 

  829.   CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopSolidNumber,CMISS_CONTROL_LOOP_NODE/), &
    830. 

  830.     & SolverSolidIndex,SolverSolid,Err)
    831. 

  831.   CALL CMISSSolver_SolverEquationsGet(SolverSolid,SolverEquationsSolid,Err)
    832. 

  832.   CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsSolid,CMISS_SOLVER_SPARSE_MATRICES,Err)
    833. 

  833.   CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsSolid,EquationsSetSolid,EquationsSetIndex,Err)
    834. 

  834.   !
    835. 

  835.   !Get the Darcy solver equations
    836. 

  836.   CALL CMISSProblem_SolverGet(Problem,(/ControlLoopSubiterationNumber,ControlLoopFluidNumber,CMISS_CONTROL_LOOP_NODE/), &
    837. 

  837.     & SolverDarcyIndex,LinearSolverDarcy,Err)
    838. 

  838.   CALL CMISSSolver_SolverEquationsGet(LinearSolverDarcy,SolverEquationsDarcy,Err)
    839. 

  839.   CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsDarcy,CMISS_SOLVER_SPARSE_MATRICES,Err)
    840. 

  840.   CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsDarcy,EquationsSetDarcy,EquationsSetIndex,Err)
    841. 

  841.   !
    842. 

  842.   CALL CMISSProblem_SolverEquationsCreateFinish(Problem,Err)
    843. 

  843. 

  844.   !
    845. 

  845.   !================================================================================================================================
    846. 

  846.   !
    847. 

  847. 

  848.   !------------------------------------
    849. 

  849.   ! ASSIGN BOUNDARY CONDITIONS - SOLID (absolute nodal parameters)
    850. 

  850.   !Solid is computed in absolute position, rather than displacement. Thus BCs for absolute position
    851. 

  851.   CALL CMISSBoundaryConditions_Initialise(BoundaryConditionsSolid,Err)
    852. 

  852.   CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsSolid,BoundaryConditionsSolid,Err)
    853. 

  853. 

  854. !   !Get surfaces
    855. 

  855. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_FRONT_SURFACE, &
    856. 

  856. !     & Face1Nodes,FaceXi(1),Err)
    857. 

  857. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BACK_SURFACE, &
    858. 

  858. !     & Face2Nodes,FaceXi(2),Err)
    859. 

  859. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_RIGHT_SURFACE, &
    860. 

  860. !     & Face3Nodes,FaceXi(3),Err)
    861. 

  861. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_LEFT_SURFACE, &
    862. 

  862. !     & Face4Nodes,FaceXi(4),Err)
    863. 

  863. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_TOP_SURFACE, &
    864. 

  864. !     & Face5Nodes,FaceXi(5),Err)
    865. 

  865. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,SolidDisplMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BOTTOM_SURFACE, &
    866. 

  866. !     & Face6Nodes,FaceXi(6),Err)
    867. 

  867. 

  868. ! ! write(*,'("Face1Nodes=[",100I3,"]")') Face1Nodes
    869. 

  869. ! ! write(*,'("Face2Nodes=[",100I3,"]")') Face2Nodes
    870. 

  870. ! ! write(*,'("Face3Nodes=[",100I3,"]")') Face3Nodes
    871. 

  871. ! ! write(*,'("Face4Nodes=[",100I3,"]")') Face4Nodes
    872. 

  872. ! ! write(*,'("Face5Nodes=[",100I3,"]")') Face5Nodes
    873. 

  873. ! ! write(*,'("Face6Nodes=[",100I3,"]")') Face6Nodes
    874. 

  874. ! ! write(*,'("FaceXi=[",100I3,"]")') FaceXi
    875. 

  875. !
    876. 

  876. allocate(Face1Nodes(21),Face2Nodes(21),Face3Nodes(21),Face4Nodes(21),Face5Nodes(9),Face6Nodes(9))
    877. 

  877. Face1Nodes=[  1, 17,  2, 22, 23, 24,  5, 31,  6, 36, 37, 38,  9, 45, 10, 50, 51, 52, 13, 59, 14]
    878. 

  878. Face2Nodes=[  3, 21,  4, 28, 29, 30,  7, 35,  8, 42, 43, 44, 11, 49, 12, 56, 57, 58, 15, 63, 16]
    879. 

  879. Face3Nodes=[  2, 20,  4, 24, 27, 30,  6, 34,  8, 38, 41, 44, 10, 48, 12, 52, 55, 58, 14, 62, 16]
    880. 

  880. Face4Nodes=[  1, 18,  3, 22, 25, 28,  5, 32,  7, 36, 39, 42,  9, 46, 11, 50, 53, 56, 13, 60, 15]
    881. 

  881. Face5Nodes=[ 13, 59, 14, 60, 61, 62, 15, 63, 16]
    882. 

  882. Face6Nodes=[  1, 17,  2, 18, 19, 20,  3, 21,  4]
    883. 

  883. FaceXi=[ -2,  2,  1, -1,  3, -3]
    884. 

  884. 

  885. 

  886.   ! Fix the bottom in z direction
    887. 

  887.   DO NN=1,SIZE(Face6Nodes,1)
    888. 

  888.     NODE=Face6Nodes(NN)
    889. 

  889.     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    890. 

  890.     IF(NodeDomain==ComputationalNodeNumber) THEN
    891. 

  891.       CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,3, &
    892. 

  892.         & ZCoord,Err)
    893. 

  893.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,3, &
    894. 

  894.         & CMISS_BOUNDARY_CONDITION_FIXED,ZCoord,Err)
    895. 

  895.       WRITE(*,*) "FIXING NODE",NODE,"AT BOTTOM IN Z DIRECTION"
    896. 

  896.     ENDIF
    897. 

  897.   ENDDO
    898. 

  898. 

  899.   ! Fix the top in z direction
    900. 

  900.   DO NN=1,SIZE(Face5Nodes,1)
    901. 

  901.     NODE=Face5Nodes(NN)
    902. 

  902.     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    903. 

  903.     IF(NodeDomain==ComputationalNodeNumber) THEN
    904. 

  904.       CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,3, &
    905. 

  905.         & ZCoord,Err)
    906. 

  906.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,3, &
    907. 

  907.         & CMISS_BOUNDARY_CONDITION_FIXED,ZCoord,Err)
    908. 

  908.       WRITE(*,*) "FIXING NODE",NODE,"AT TOP IN Z DIRECTION"
    909. 

  909.     ENDIF
    910. 

  910.   ENDDO
    911. 

  911. 

  912.   !Fix more nodes at the bottom to stop free body motion
    913. 

  913.   X_FIXED=.FALSE.
    914. 

  914.   Y_FIXED=.FALSE.
    915. 

  915.   DO NN=1,SIZE(Face6Nodes,1)
    916. 

  916.     NODE=Face6Nodes(NN)
    917. 

  917.     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    918. 

  918.     IF(NodeDomain==ComputationalNodeNumber) THEN
    919. 

  919.       CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,1, &
    920. 

  920.         & XCoord,Err)
    921. 

  921.       CALL CMISSField_ParameterSetGetNode(GeometricFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,1,NODE,2, &
    922. 

  922.         & YCoord,Err)
    923. 

  923. 

  924.       !Fix Origin displacement in x and y (z already fixed)
    925. 

  925.       IF(ABS(XCoord)<1.0E-6_CMISSDP.AND.ABS(YCoord)<1.0E-6_CMISSDP) THEN
    926. 

  926.         CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,1, &
    927. 

  927.           & CMISS_BOUNDARY_CONDITION_FIXED,XCoord,Err)
    928. 

  928.         CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,2, &
    929. 

  929.           & CMISS_BOUNDARY_CONDITION_FIXED,YCoord,Err)
    930. 

  930.         WRITE(*,*) "FIXING ORIGIN NODE",NODE,"IN X AND Y DIRECTION"
    931. 

  931.         X_FIXED=.TRUE.
    932. 

  932.         Y_FIXED=.TRUE.
    933. 

  933.       ENDIF
    934. 

  934. 

  935.       !Fix nodal displacements at (X_DIM,0) in y
    936. 

  936.       IF(ABS(XCoord - X_DIM)<1.0E-6_CMISSDP .AND. ABS(YCoord)<1.0E-6_CMISSDP) THEN
    937. 

  937.         CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,2, &
    938. 

  938.           & CMISS_BOUNDARY_CONDITION_FIXED,YCoord,Err)
    939. 

  939.         WRITE(*,*) "FIXING NODES",NODE,"AT (X_DIM,0) IN Y DIRECTION"
    940. 

  940.         Y_FIXED=.TRUE.
    941. 

  941.       ENDIF
    942. 

  942. 

  943.       !Fix nodal displacements at (0,Y_DIM) in x
    944. 

  944.       IF(ABS(XCoord)<1.0E-6_CMISSDP .AND. ABS(YCoord - Y_DIM)<1.0E-6_CMISSDP) THEN
    945. 

  945.         CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsSolid,DependentFieldSolid,CMISS_FIELD_U_VARIABLE_TYPE,1,1,NODE,1, &
    946. 

  946.           & CMISS_BOUNDARY_CONDITION_FIXED,XCoord,Err)
    947. 

  947.         WRITE(*,*) "FIXING NODES",NODE,"AT (0,Y_DIM) IN X DIRECTION"
    948. 

  948.         X_FIXED=.TRUE.
    949. 

  949.       ENDIF
    950. 

  950. 

  951.     ENDIF
    952. 

  952.   ENDDO
    953. 

  953. !   CALL MPI_REDUCE(X_FIXED,X_OKAY,1,MPI_LOGICAL,MPI_LOR,0,MPI_COMM_WORLD,MPI_IERROR)
    954. 

  954. !   CALL MPI_REDUCE(Y_FIXED,Y_OKAY,1,MPI_LOGICAL,MPI_LOR,0,MPI_COMM_WORLD,MPI_IERROR)
    955. 

  955. !   IF(ComputationalNodeNumber==0) THEN
    956. 

  956. !     IF(.NOT.(X_OKAY.AND.Y_OKAY)) THEN
    957. 

  957. !       WRITE(*,*) "Free body motion could not be prevented!"
    958. 

  958. !       CALL CMISSFinalise(Err)
    959. 

  959. !       STOP
    960. 

  960. !     ENDIF
    961. 

  961. !   ENDIF
    962. 

  962. 

  963.   CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquationsSolid,Err)
    964. 

  964.   !------------------------------------
    965. 

  965. 

  966. 

  967.   !------------------------------------
    968. 

  968.   ! ASSIGN BOUNDARY CONDITIONS - FLUID
    969. 

  969.   CALL CMISSBoundaryConditions_Initialise(BoundaryConditionsDarcy,Err)
    970. 

  970.   CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsDarcy,BoundaryConditionsDarcy,Err)
    971. 

  971. 

  972. !   !Get surfaces
    973. 

  973. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_FRONT_SURFACE, &
    974. 

  974. !     & Face7Nodes,FaceXi(1),Err)
    975. 

  975. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BACK_SURFACE, &
    976. 

  976. !     & Face8Nodes,FaceXi(2),Err)
    977. 

  977. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_RIGHT_SURFACE, &
    978. 

  978. !     & Face9Nodes,FaceXi(3),Err)
    979. 

  979. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_LEFT_SURFACE, &
    980. 

  980. !     & Face10Nodes,FaceXi(4),Err)
    981. 

  981. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_TOP_SURFACE, &
    982. 

  982. !     & Face11Nodes,FaceXi(5),Err)
    983. 

  983. !   CALL CMISSGeneratedMesh_SurfaceGet(GeneratedMesh,DarcyVelMeshComponentNumber,CMISS_GENERATED_MESH_REGULAR_BOTTOM_SURFACE, &
    984. 

  984. !     & Face12Nodes,FaceXi(6),Err)
    985. 

  985. 

  986. ! ! write(*,'("Face7Nodes=[",100I3,"]")') Face7Nodes
    987. 

  987. ! ! write(*,'("Face8Nodes=[",100I3,"]")') Face8Nodes
    988. 

  988. ! ! write(*,'("Face9Nodes=[",100I3,"]")') Face9Nodes
    989. 

  989. ! ! write(*,'("Face10Nodes=[",100I3,"]")') Face10Nodes
    990. 

  990. ! ! write(*,'("Face11Nodes=[",100I3,"]")') Face11Nodes
    991. 

  991. ! ! write(*,'("Face12Nodes=[",100I3,"]")') Face12Nodes
    992. 

  992. ! ! write(*,'("FaceXi=[",100I3,"]")') FaceXi
    993. 

  993. 

  994. allocate(Face7Nodes(8),Face8Nodes(8),Face9Nodes(8),Face10Nodes(8),Face11Nodes(4),Face12Nodes(4))
    995. 

  995. Face7Nodes=[  1,  2,  5,  6,  9, 10, 13, 14]
    996. 

  996. Face8Nodes=[  3,  4,  7,  8, 11, 12, 15, 16]
    997. 

  997. Face9Nodes=[  2,  4,  6,  8, 10, 12, 14, 16]
    998. 

  998. Face10Nodes=[  1,  3,  5,  7,  9, 11, 13, 15]
    999. 

  999. Face11Nodes=[ 13, 14, 15, 16]
    1000. 

  1000. Face12Nodes=[  1,  2,  3,  4]
    1001. 

  1001. FaceXi=[ -2,  2,  1, -1,  3, -3]
    1002. 

  1002. 

  1003. 

  1004.   ! At the top impose Darcy velocity in z direction
    1005. 

  1005.   DO NN=1,SIZE(Face11Nodes,1)
    1006. 

  1006.     NODE=Face11Nodes(NN)
    1007. 

  1007. !     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    1008. 

  1008. !     IF(NodeDomain==ComputationalNodeNumber) THEN
    1009. 

  1009.       VALUE = -2.0_CMISSDP
    1010. 

  1010.       COMPONENT_NUMBER = 3
    1011. 

  1011.       write(*,*)'Marker 0'
    1012. 

  1012.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsDarcy,DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,1,NODE, &
    1013. 

  1013.         & COMPONENT_NUMBER, &
    1014. 

  1014.         & CMISS_BOUNDARY_CONDITION_FIXED,VALUE,Err)
    1015. 

  1015.       WRITE(*,*) "SPECIFIED INFLOW AT NODE",NODE,"IN Z DIRECTION"
    1016. 

  1016. 

  1017. !       CALL CMISSField_ParameterSetGetNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,NODE,1,XCoord,Err)
    1018. 

  1018. !       CALL CMISSField_ParameterSetGetNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,NODE,2,YCoord,Err)
    1019. 

  1019. !       CALL CMISSField_ParameterSetGetNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,1,NODE,3,ZCoord,Err)
    1020. 

  1020. !       WRITE(*,*) "XCoord, YCoord, ZCoord = ",XCoord, YCoord, ZCoord
    1021. 

  1021. !     ENDIF
    1022. 

  1022.   ENDDO
    1023. 

  1023. 

  1024.   !All other faces are impermeable
    1025. 

  1025.   DO NN=1,SIZE(Face7Nodes,1)
    1026. 

  1026.     NODE=Face7Nodes(NN)
    1027. 

  1027. !     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    1028. 

  1028. !     IF(NodeDomain==ComputationalNodeNumber) THEN
    1029. 

  1029.       VALUE = 0.0_CMISSDP
    1030. 

  1030.       COMPONENT_NUMBER = 1
    1031. 

  1031.       write(*,*)'Marker 1'
    1032. 

  1032.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsDarcy,DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,1,NODE, &
    1033. 

  1033.         & COMPONENT_NUMBER, &
    1034. 

  1034.         & CMISS_BOUNDARY_CONDITION_FIXED,VALUE,Err)
    1035. 

  1035.       WRITE(*,*) "SPECIFIED IMPERMEABLE WALL AT NODE",NODE,"IN X DIRECTION"
    1036. 

  1036. !     ENDIF
    1037. 

  1037.   ENDDO
    1038. 

  1038. 

  1039.   DO NN=1,SIZE(Face8Nodes,1)
    1040. 

  1040.     NODE=Face8Nodes(NN)
    1041. 

  1041. !     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    1042. 

  1042. !     IF(NodeDomain==ComputationalNodeNumber) THEN
    1043. 

  1043.       VALUE = 0.0_CMISSDP
    1044. 

  1044.       COMPONENT_NUMBER = 1
    1045. 

  1045.       write(*,*)'Marker 2'
    1046. 

  1046.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsDarcy,DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,1,NODE, &
    1047. 

  1047.         & COMPONENT_NUMBER, &
    1048. 

  1048.         & CMISS_BOUNDARY_CONDITION_FIXED,VALUE,Err)
    1049. 

  1049.       WRITE(*,*) "SPECIFIED IMPERMEABLE WALL AT NODE",NODE,"IN X DIRECTION"
    1050. 

  1050. !     ENDIF
    1051. 

  1051.   ENDDO
    1052. 

  1052. 

  1053.   DO NN=1,SIZE(Face9Nodes,1)
    1054. 

  1054.     NODE=Face9Nodes(NN)
    1055. 

  1055. !     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    1056. 

  1056. !     IF(NodeDomain==ComputationalNodeNumber) THEN
    1057. 

  1057.       VALUE = 0.0_CMISSDP
    1058. 

  1058.       COMPONENT_NUMBER = 2
    1059. 

  1059.       write(*,*)'Marker 3'
    1060. 

  1060.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsDarcy,DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,1,NODE, &
    1061. 

  1061.         & COMPONENT_NUMBER, &
    1062. 

  1062.         & CMISS_BOUNDARY_CONDITION_FIXED,VALUE,Err)
    1063. 

  1063.       WRITE(*,*) "SPECIFIED IMPERMEABLE WALL AT NODE",NODE,"IN Y DIRECTION"
    1064. 

  1064. !     ENDIF
    1065. 

  1065.   ENDDO
    1066. 

  1066. 

  1067.   DO NN=1,SIZE(Face10Nodes,1)
    1068. 

  1068.     NODE=Face10Nodes(NN)
    1069. 

  1069. !     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    1070. 

  1070. !     IF(NodeDomain==ComputationalNodeNumber) THEN
    1071. 

  1071.       VALUE = 0.0_CMISSDP
    1072. 

  1072.       COMPONENT_NUMBER = 2
    1073. 

  1073.       write(*,*)'Marker 4'
    1074. 

  1074.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsDarcy,DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,1,NODE, &
    1075. 

  1075.         & COMPONENT_NUMBER, &
    1076. 

  1076.         & CMISS_BOUNDARY_CONDITION_FIXED,VALUE,Err)
    1077. 

  1077.       WRITE(*,*) "SPECIFIED IMPERMEABLE WALL AT NODE",NODE,"IN Y DIRECTION"
    1078. 

  1078. !     ENDIF
    1079. 

  1079.   ENDDO
    1080. 

  1080. 

  1081.   DO NN=1,SIZE(Face12Nodes,1)
    1082. 

  1082.     NODE=Face12Nodes(NN)
    1083. 

  1083. !     CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE,1,NodeDomain,Err)
    1084. 

  1084. !     IF(NodeDomain==ComputationalNodeNumber) THEN
    1085. 

  1085.       VALUE = 0.0_CMISSDP
    1086. 

  1086.       COMPONENT_NUMBER = 3
    1087. 

  1087.       write(*,*)'Marker 5'
    1088. 

  1088.       CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsDarcy,DependentFieldSolid,CMISS_FIELD_V_VARIABLE_TYPE,1,1,NODE, &
    1089. 

  1089.         & COMPONENT_NUMBER, &
    1090. 

  1090.         & CMISS_BOUNDARY_CONDITION_FIXED,VALUE,Err)
    1091. 

  1091.       WRITE(*,*) "SPECIFIED IMPERMEABLE WALL AT NODE",NODE,"IN Z DIRECTION"
    1092. 

  1092. !     ENDIF
    1093. 

  1093.   ENDDO
    1094. 

  1094. 

  1095.   CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquationsDarcy,Err)
    1096. 

  1096. 

  1097.   !
    1098. 

  1098.   !================================================================================================================================
    1099. 

  1099.   !
    1100. 

  1100. 

  1101.   !RUN SOLVERS
    1102. 

  1102. 

  1103.   !Turn of PETSc error handling
    1104. 

  1104.   !CALL PETSC_ERRORHANDLING_SET_ON(ERR,ERROR,*999)
    1105. 

  1105. 

  1106.   !Solve the problem
    1107. 

  1107.   WRITE(*,'(A)') "Solving problem..."
    1108. 

  1108.   CALL CMISSProblem_Solve(Problem,Err)
    1109. 

  1109.   WRITE(*,'(A)') "Problem solved!"
    1110. 

  1110. 

  1111. 

  1112.   !
    1113. 

  1113.   !================================================================================================================================
    1114. 

  1114.   !
    1115. 

  1115. 

  1116.   !OUTPUT
    1117. 

  1117. 

  1118.   EXPORT_FIELD_IO=.FALSE.
    1119. 

  1119.   IF(EXPORT_FIELD_IO) THEN
    1120. 

  1120.     WRITE(*,'(A)') "Exporting fields..."
    1121. 

  1121.     CALL CMISSFields_Initialise(Fields,Err)
    1122. 

  1122.     CALL CMISSFields_Create(Region,Fields,Err)
    1123. 

  1123.     CALL CMISSFields_NodesExport(Fields,"FiniteElasticityDarcy","FORTRAN",Err)
    1124. 

  1124.     CALL CMISSFields_ElementsExport(Fields,"FiniteElasticityDarcy","FORTRAN",Err)
    1125. 

  1125.     CALL CMISSFields_Finalise(Fields,Err)
    1126. 

  1126.     WRITE(*,'(A)') "Field exported!"
    1127. 

  1127.   ENDIF
    1128. 

  1128. 

  1129. 

  1130.   !Finialise CMISS
    1131. 

  1131. !   CALL CMISSFinalise(Err)
    1132. 

  1132. 

  1133.   WRITE(*,'(A)') "Program successfully completed."
    1134. 

  1134. 

  1135.   STOP
    1136. 

  1136. 

  1137. END PROGRAM FINITEELASTICITYDARCYIOEXAMPLE
    1138. 

  1138.