/FluidMechanics/Stokes/RoutineCheck/Dynamic/src/DynamicExample.f90

http://github.com/xyan075/examples · Fortran Modern · 800 lines · 463 code · 104 blank · 233 comment · 0 complexity · ce21518a693177ae087e6cee65794220 MD5 · raw file

  1. !> \file
  2. !> \author Sebastian Krittian
  3. !> \brief This is an example program to solve a dynamic Stokes 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. !> \example FluidMechanics/Stokes/RoutineCheck/Dynamic/src/DynamicExample.f90
  42. !! Example program to solve a dynamic Stokes equation using OpenCMISS calls.
  43. !! \htmlinclude FluidMechanics/Stokes/RoutineCheck/Dynamic/history.html
  44. !!
  45. !<
  46. !> Main program
  47. PROGRAM STOKESDYNAMICEXAMPLE
  48. !
  49. !================================================================================================================================
  50. !
  51. !PROGRAM LIBRARIES
  52. USE OPENCMISS
  53. USE FLUID_MECHANICS_IO_ROUTINES
  54. USE MPI
  55. #ifdef WIN32
  56. USE IFQWINCMISS
  57. #endif
  58. !
  59. !================================================================================================================================
  60. !
  61. !PROGRAM VARIABLES AND TYPES
  62. IMPLICIT NONE
  63. INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldUserNumber=1337
  64. TYPE(CMISSFieldType) :: EquationsSetField
  65. !Test program parameters
  66. INTEGER(CMISSIntg), PARAMETER :: CoordinateSystemUserNumber=1
  67. INTEGER(CMISSIntg), PARAMETER :: RegionUserNumber=2
  68. INTEGER(CMISSIntg), PARAMETER :: MeshUserNumber=3
  69. INTEGER(CMISSIntg), PARAMETER :: DecompositionUserNumber=4
  70. INTEGER(CMISSIntg), PARAMETER :: GeometricFieldUserNumber=5
  71. INTEGER(CMISSIntg), PARAMETER :: DependentFieldUserNumberStokes=6
  72. INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberStokes=7
  73. INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumberStokes=8
  74. INTEGER(CMISSIntg), PARAMETER :: ProblemUserNumber=9
  75. INTEGER(CMISSIntg), PARAMETER :: DomainUserNumber=2
  76. INTEGER(CMISSIntg), PARAMETER :: SolverStokesUserNumber=1
  77. INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberStokesMu=1
  78. INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberStokesRho=2
  79. !Program types
  80. TYPE(EXPORT_CONTAINER):: CM
  81. !Program variables
  82. INTEGER(CMISSIntg) :: NUMBER_OF_DIMENSIONS
  83. INTEGER(CMISSIntg) :: BASIS_TYPE
  84. INTEGER(CMISSIntg) :: BASIS_NUMBER_SPACE
  85. INTEGER(CMISSIntg) :: BASIS_NUMBER_VELOCITY
  86. INTEGER(CMISSIntg) :: BASIS_NUMBER_PRESSURE
  87. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_SPACE
  88. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_VELOCITY
  89. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_PRESSURE
  90. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_SPACE
  91. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_VELOCITY
  92. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_PRESSURE
  93. INTEGER(CMISSIntg) :: MESH_NUMBER_OF_COMPONENTS
  94. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_SPACE
  95. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_VELOCITY
  96. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_PRESSURE
  97. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_SPACE
  98. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_VELOCITY
  99. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_PRESSURE
  100. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_SPACE
  101. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_VELOCITY
  102. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_PRESSURE
  103. INTEGER(CMISSIntg) :: TOTAL_NUMBER_OF_NODES
  104. INTEGER(CMISSIntg) :: TOTAL_NUMBER_OF_ELEMENTS
  105. INTEGER(CMISSIntg) :: MAXIMUM_ITERATIONS
  106. INTEGER(CMISSIntg) :: RESTART_VALUE
  107. ! INTEGER(CMISSIntg) :: MPI_IERROR
  108. INTEGER(CMISSIntg) :: NUMBER_OF_FIXED_WALL_NODES_STOKES
  109. INTEGER(CMISSIntg) :: NUMBER_OF_INLET_WALL_NODES_STOKES
  110. INTEGER(CMISSIntg) :: EQUATIONS_STOKES_OUTPUT
  111. INTEGER(CMISSIntg) :: COMPONENT_NUMBER
  112. INTEGER(CMISSIntg) :: NODE_NUMBER
  113. INTEGER(CMISSIntg) :: ELEMENT_NUMBER
  114. INTEGER(CMISSIntg) :: NODE_COUNTER
  115. INTEGER(CMISSIntg) :: CONDITION
  116. INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_STOKES_OUTPUT_FREQUENCY
  117. INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_STOKES_OUTPUT_TYPE
  118. INTEGER(CMISSIntg) :: LINEAR_SOLVER_STOKES_OUTPUT_TYPE
  119. INTEGER, ALLOCATABLE, DIMENSION(:):: FIXED_WALL_NODES_STOKES
  120. INTEGER, ALLOCATABLE, DIMENSION(:):: INLET_WALL_NODES_STOKES
  121. REAL(CMISSDP) :: INITIAL_FIELD_STOKES(3)
  122. REAL(CMISSDP) :: BOUNDARY_CONDITIONS_STOKES(3)
  123. REAL(CMISSDP) :: DIVERGENCE_TOLERANCE
  124. REAL(CMISSDP) :: RELATIVE_TOLERANCE
  125. REAL(CMISSDP) :: ABSOLUTE_TOLERANCE
  126. REAL(CMISSDP) :: LINESEARCH_ALPHA
  127. REAL(CMISSDP) :: VALUE
  128. REAL(CMISSDP) :: MU_PARAM_STOKES
  129. REAL(CMISSDP) :: RHO_PARAM_STOKES
  130. REAL(CMISSDP) :: DYNAMIC_SOLVER_STOKES_START_TIME
  131. REAL(CMISSDP) :: DYNAMIC_SOLVER_STOKES_STOP_TIME
  132. REAL(CMISSDP) :: DYNAMIC_SOLVER_STOKES_THETA
  133. REAL(CMISSDP) :: DYNAMIC_SOLVER_STOKES_TIME_INCREMENT
  134. LOGICAL :: EXPORT_FIELD_IO
  135. LOGICAL :: LINEAR_SOLVER_STOKES_DIRECT_FLAG
  136. LOGICAL :: FIXED_WALL_NODES_STOKES_FLAG
  137. LOGICAL :: INLET_WALL_NODES_STOKES_FLAG
  138. !CMISS variables
  139. !Regions
  140. TYPE(CMISSRegionType) :: Region
  141. TYPE(CMISSRegionType) :: WorldRegion
  142. !Coordinate systems
  143. TYPE(CMISSCoordinateSystemType) :: CoordinateSystem
  144. TYPE(CMISSCoordinateSystemType) :: WorldCoordinateSystem
  145. !Basis
  146. TYPE(CMISSBasisType) :: BasisSpace
  147. TYPE(CMISSBasisType) :: BasisVelocity
  148. TYPE(CMISSBasisType) :: BasisPressure
  149. !Nodes
  150. TYPE(CMISSNodesType) :: Nodes
  151. !Elements
  152. TYPE(CMISSMeshElementsType) :: MeshElementsSpace
  153. TYPE(CMISSMeshElementsType) :: MeshElementsVelocity
  154. TYPE(CMISSMeshElementsType) :: MeshElementsPressure
  155. !Meshes
  156. TYPE(CMISSMeshType) :: Mesh
  157. !Decompositions
  158. TYPE(CMISSDecompositionType) :: Decomposition
  159. !Fields
  160. TYPE(CMISSFieldsType) :: Fields
  161. !Field types
  162. TYPE(CMISSFieldType) :: GeometricField
  163. TYPE(CMISSFieldType) :: DependentFieldStokes
  164. TYPE(CMISSFieldType) :: MaterialsFieldStokes
  165. !Boundary conditions
  166. TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsStokes
  167. !Equations sets
  168. TYPE(CMISSEquationsSetType) :: EquationsSetStokes
  169. !Equations
  170. TYPE(CMISSEquationsType) :: EquationsStokes
  171. !Problems
  172. TYPE(CMISSProblemType) :: Problem
  173. !Control loops
  174. TYPE(CMISSControlLoopType) :: ControlLoop
  175. !Solvers
  176. TYPE(CMISSSolverType) :: DynamicSolverStokes
  177. TYPE(CMISSSolverType) :: LinearSolverStokes
  178. !Solver equations
  179. TYPE(CMISSSolverEquationsType) :: SolverEquationsStokes
  180. #ifdef WIN32
  181. !Quickwin type
  182. LOGICAL :: QUICKWIN_STATUS=.FALSE.
  183. TYPE(WINDOWCONFIG) :: QUICKWIN_WINDOW_CONFIG
  184. #endif
  185. !Generic CMISS variables
  186. INTEGER(CMISSIntg) :: NumberOfComputationalNodes,ComputationalNodeNumber,BoundaryNodeDomain
  187. INTEGER(CMISSIntg) :: EquationsSetIndex
  188. INTEGER(CMISSIntg) :: Err
  189. #ifdef WIN32
  190. !Initialise QuickWin
  191. QUICKWIN_WINDOW_CONFIG%TITLE="General Output" !Window title
  192. QUICKWIN_WINDOW_CONFIG%NUMTEXTROWS=-1 !Max possible number of rows
  193. QUICKWIN_WINDOW_CONFIG%MODE=QWIN$SCROLLDOWN
  194. !Set the window parameters
  195. QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
  196. !If attempt fails set with system estimated values
  197. IF(.NOT.QUICKWIN_STATUS) QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
  198. #endif
  199. !
  200. !================================================================================================================================
  201. !
  202. !INITIALISE OPENCMISS
  203. CALL CMISSInitialise(WorldCoordinateSystem,WorldRegion,Err)
  204. CALL CMISSErrorHandlingModeSet(CMISS_ERRORS_TRAP_ERROR,Err)
  205. !
  206. !================================================================================================================================
  207. !
  208. !CHECK COMPUTATIONAL NODE
  209. !Get the computational nodes information
  210. CALL CMISSComputationalNumberOfNodesGet(NumberOfComputationalNodes,Err)
  211. CALL CMISSComputationalNodeNumberGet(ComputationalNodeNumber,Err)
  212. !
  213. !================================================================================================================================
  214. !
  215. !PROBLEM CONTROL PANEL
  216. !Import cmHeart mesh information
  217. CALL FLUID_MECHANICS_IO_READ_CMHEART(CM,Err)
  218. BASIS_NUMBER_SPACE=CM%ID_M
  219. BASIS_NUMBER_VELOCITY=CM%ID_V
  220. BASIS_NUMBER_PRESSURE=CM%ID_P
  221. NUMBER_OF_DIMENSIONS=CM%D
  222. BASIS_TYPE=CM%IT_T
  223. BASIS_XI_INTERPOLATION_SPACE=CM%IT_M
  224. BASIS_XI_INTERPOLATION_VELOCITY=CM%IT_V
  225. BASIS_XI_INTERPOLATION_PRESSURE=CM%IT_P
  226. NUMBER_OF_NODES_SPACE=CM%N_M
  227. NUMBER_OF_NODES_VELOCITY=CM%N_V
  228. NUMBER_OF_NODES_PRESSURE=CM%N_P
  229. TOTAL_NUMBER_OF_NODES=CM%N_T
  230. TOTAL_NUMBER_OF_ELEMENTS=CM%E_T
  231. NUMBER_OF_ELEMENT_NODES_SPACE=CM%EN_M
  232. NUMBER_OF_ELEMENT_NODES_VELOCITY=CM%EN_V
  233. NUMBER_OF_ELEMENT_NODES_PRESSURE=CM%EN_P
  234. !Set initial values
  235. INITIAL_FIELD_STOKES(1)=0.0_CMISSDP
  236. INITIAL_FIELD_STOKES(2)=0.0_CMISSDP
  237. INITIAL_FIELD_STOKES(3)=0.0_CMISSDP
  238. !Set boundary conditions
  239. FIXED_WALL_NODES_STOKES_FLAG=.TRUE.
  240. INLET_WALL_NODES_STOKES_FLAG=.TRUE.
  241. IF(FIXED_WALL_NODES_STOKES_FLAG) THEN
  242. NUMBER_OF_FIXED_WALL_NODES_STOKES=80
  243. ALLOCATE(FIXED_WALL_NODES_STOKES(NUMBER_OF_FIXED_WALL_NODES_STOKES))
  244. FIXED_WALL_NODES_STOKES=(/1,2,3,4,5,7,9,10,11,12,13,14,17,20,24,28,29,30,31,32,33,34,35,37,39, &
  245. & 41,44,46,47,48,50,51,52,53,54,57,60,64,65,66,67,68,70,72,74,76,77,78,79,80,83,86, &
  246. & 89,90,91,92,93,94,95,97,99,101,102,103,104,105,106,107,108,111,114,115,116,117,118, &
  247. & 120,122,123,124,125/)
  248. ENDIF
  249. IF(INLET_WALL_NODES_STOKES_FLAG) THEN
  250. NUMBER_OF_INLET_WALL_NODES_STOKES=9
  251. ALLOCATE(INLET_WALL_NODES_STOKES(NUMBER_OF_INLET_WALL_NODES_STOKES))
  252. INLET_WALL_NODES_STOKES=(/6,15,16,23,36,42,81,82,96/)
  253. !Set initial boundary conditions
  254. BOUNDARY_CONDITIONS_STOKES(1)=0.0_CMISSDP
  255. BOUNDARY_CONDITIONS_STOKES(2)=1.0_CMISSDP
  256. BOUNDARY_CONDITIONS_STOKES(3)=0.0_CMISSDP
  257. ENDIF
  258. !Set material parameters
  259. MU_PARAM_STOKES=1.0_CMISSDP
  260. RHO_PARAM_STOKES=1.0_CMISSDP
  261. !Set interpolation parameters
  262. BASIS_XI_GAUSS_SPACE=3
  263. BASIS_XI_GAUSS_VELOCITY=3
  264. BASIS_XI_GAUSS_PRESSURE=3
  265. !Set output parameter
  266. !(NoOutput/ProgressOutput/TimingOutput/SolverOutput/SolverMatrixOutput)
  267. DYNAMIC_SOLVER_STOKES_OUTPUT_TYPE=CMISS_SOLVER_NO_OUTPUT
  268. LINEAR_SOLVER_STOKES_OUTPUT_TYPE=CMISS_SOLVER_NO_OUTPUT
  269. !(NoOutput/TimingOutput/MatrixOutput/ElementOutput)
  270. EQUATIONS_STOKES_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT
  271. !Set time parameter
  272. DYNAMIC_SOLVER_STOKES_START_TIME=0.0_CMISSDP
  273. DYNAMIC_SOLVER_STOKES_STOP_TIME=5.0_CMISSDP
  274. DYNAMIC_SOLVER_STOKES_TIME_INCREMENT=1.0_CMISSDP
  275. DYNAMIC_SOLVER_STOKES_THETA=1.0_CMISSDP
  276. !Set result output parameter
  277. DYNAMIC_SOLVER_STOKES_OUTPUT_FREQUENCY=1
  278. !Set solver parameters
  279. LINEAR_SOLVER_STOKES_DIRECT_FLAG=.FALSE.
  280. RELATIVE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-05_CMISSDP
  281. ABSOLUTE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-10_CMISSDP
  282. DIVERGENCE_TOLERANCE=1.0E20 !default: 1.0E5
  283. MAXIMUM_ITERATIONS=100000 !default: 100000
  284. RESTART_VALUE=3000 !default: 30
  285. LINESEARCH_ALPHA=1.0
  286. !
  287. !================================================================================================================================
  288. !
  289. !COORDINATE SYSTEM
  290. !Start the creation of a new RC coordinate system
  291. CALL CMISSCoordinateSystem_Initialise(CoordinateSystem,Err)
  292. CALL CMISSCoordinateSystem_CreateStart(CoordinateSystemUserNumber,CoordinateSystem,Err)
  293. !Set the coordinate system dimension
  294. CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,NUMBER_OF_DIMENSIONS,Err)
  295. !Finish the creation of the coordinate system
  296. CALL CMISSCoordinateSystem_CreateFinish(CoordinateSystem,Err)
  297. !
  298. !================================================================================================================================
  299. !
  300. !REGION
  301. !Start the creation of a new region
  302. CALL CMISSRegion_Initialise(Region,Err)
  303. CALL CMISSRegion_CreateStart(RegionUserNumber,WorldRegion,Region,Err)
  304. !Set the regions coordinate system as defined above
  305. CALL CMISSRegion_CoordinateSystemSet(Region,CoordinateSystem,Err)
  306. !Finish the creation of the region
  307. CALL CMISSRegion_CreateFinish(Region,Err)
  308. !
  309. !================================================================================================================================
  310. !
  311. !BASES
  312. !Start the creation of new bases
  313. MESH_NUMBER_OF_COMPONENTS=1
  314. CALL CMISSBasis_Initialise(BasisSpace,Err)
  315. CALL CMISSBasis_CreateStart(BASIS_NUMBER_SPACE,BasisSpace,Err)
  316. !Set the basis type (Lagrange/Simplex)
  317. CALL CMISSBasis_TypeSet(BasisSpace,BASIS_TYPE,Err)
  318. !Set the basis xi number
  319. CALL CMISSBasis_NumberOfXiSet(BasisSpace,NUMBER_OF_DIMENSIONS,Err)
  320. !Set the basis xi interpolation and number of Gauss points
  321. IF(NUMBER_OF_DIMENSIONS==2) THEN
  322. CALL CMISSBasis_InterpolationXiSet(BasisSpace,(/BASIS_XI_INTERPOLATION_SPACE,BASIS_XI_INTERPOLATION_SPACE/),Err)
  323. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisSpace,(/BASIS_XI_GAUSS_SPACE,BASIS_XI_GAUSS_SPACE/),Err)
  324. ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
  325. CALL CMISSBasis_InterpolationXiSet(BasisSpace,(/BASIS_XI_INTERPOLATION_SPACE,BASIS_XI_INTERPOLATION_SPACE, &
  326. & BASIS_XI_INTERPOLATION_SPACE/),Err)
  327. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisSpace,(/BASIS_XI_GAUSS_SPACE,BASIS_XI_GAUSS_SPACE,BASIS_XI_GAUSS_SPACE/),Err)
  328. ENDIF
  329. !Finish the creation of the basis
  330. CALL CMISSBasis_CreateFinish(BasisSpace,Err)
  331. !Start the creation of another basis
  332. IF(BASIS_XI_INTERPOLATION_VELOCITY==BASIS_XI_INTERPOLATION_SPACE) THEN
  333. BasisVelocity=BasisSpace
  334. ELSE
  335. MESH_NUMBER_OF_COMPONENTS=MESH_NUMBER_OF_COMPONENTS+1
  336. !Initialise a new velocity basis
  337. CALL CMISSBasis_Initialise(BasisVelocity,Err)
  338. !Start the creation of a basis
  339. CALL CMISSBasis_CreateStart(BASIS_NUMBER_VELOCITY,BasisVelocity,Err)
  340. !Set the basis type (Lagrange/Simplex)
  341. CALL CMISSBasis_TypeSet(BasisVelocity,BASIS_TYPE,Err)
  342. !Set the basis xi number
  343. CALL CMISSBasis_NumberOfXiSet(BasisVelocity,NUMBER_OF_DIMENSIONS,Err)
  344. !Set the basis xi interpolation and number of Gauss points
  345. IF(NUMBER_OF_DIMENSIONS==2) THEN
  346. CALL CMISSBasis_InterpolationXiSet(BasisVelocity,(/BASIS_XI_INTERPOLATION_VELOCITY,BASIS_XI_INTERPOLATION_VELOCITY/),Err)
  347. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisVelocity,(/BASIS_XI_GAUSS_VELOCITY,BASIS_XI_GAUSS_VELOCITY/),Err)
  348. ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
  349. CALL CMISSBasis_InterpolationXiSet(BasisVelocity,(/BASIS_XI_INTERPOLATION_VELOCITY,BASIS_XI_INTERPOLATION_VELOCITY, &
  350. & BASIS_XI_INTERPOLATION_VELOCITY/),Err)
  351. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisVelocity,(/BASIS_XI_GAUSS_VELOCITY,BASIS_XI_GAUSS_VELOCITY, &
  352. & BASIS_XI_GAUSS_VELOCITY/),Err)
  353. ENDIF
  354. !Finish the creation of the basis
  355. CALL CMISSBasis_CreateFinish(BasisVelocity,Err)
  356. ENDIF
  357. !Start the creation of another basis
  358. IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_SPACE) THEN
  359. BasisPressure=BasisSpace
  360. ELSE IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_VELOCITY) THEN
  361. BasisPressure=BasisVelocity
  362. ELSE
  363. MESH_NUMBER_OF_COMPONENTS=MESH_NUMBER_OF_COMPONENTS+1
  364. !Initialise a new pressure basis
  365. CALL CMISSBasis_Initialise(BasisPressure,Err)
  366. !Start the creation of a basis
  367. CALL CMISSBasis_CreateStart(BASIS_NUMBER_PRESSURE,BasisPressure,Err)
  368. !Set the basis type (Lagrange/Simplex)
  369. CALL CMISSBasis_TypeSet(BasisPressure,BASIS_TYPE,Err)
  370. !Set the basis xi number
  371. CALL CMISSBasis_NumberOfXiSet(BasisPressure,NUMBER_OF_DIMENSIONS,Err)
  372. !Set the basis xi interpolation and number of Gauss points
  373. IF(NUMBER_OF_DIMENSIONS==2) THEN
  374. CALL CMISSBasis_InterpolationXiSet(BasisPressure,(/BASIS_XI_INTERPOLATION_PRESSURE,BASIS_XI_INTERPOLATION_PRESSURE/),Err)
  375. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisPressure,(/BASIS_XI_GAUSS_PRESSURE,BASIS_XI_GAUSS_PRESSURE/),Err)
  376. ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
  377. CALL CMISSBasis_InterpolationXiSet(BasisPressure,(/BASIS_XI_INTERPOLATION_PRESSURE,BASIS_XI_INTERPOLATION_PRESSURE, &
  378. & BASIS_XI_INTERPOLATION_PRESSURE/),Err)
  379. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisPressure,(/BASIS_XI_GAUSS_PRESSURE,BASIS_XI_GAUSS_PRESSURE, &
  380. & BASIS_XI_GAUSS_PRESSURE/),Err)
  381. ENDIF
  382. !Finish the creation of the basis
  383. CALL CMISSBasis_CreateFinish(BasisPressure,Err)
  384. ENDIF
  385. !
  386. !================================================================================================================================
  387. !
  388. !MESH
  389. !Start the creation of mesh nodes
  390. CALL CMISSNodes_Initialise(Nodes,Err)
  391. CALL CMISSMesh_Initialise(Mesh,Err)
  392. CALL CMISSNodes_CreateStart(Region,TOTAL_NUMBER_OF_NODES,Nodes,Err)
  393. CALL CMISSNodes_CreateFinish(Nodes,Err)
  394. !Start the creation of the mesh
  395. CALL CMISSMesh_CreateStart(MeshUserNumber,Region,NUMBER_OF_DIMENSIONS,Mesh,Err)
  396. !Set number of mesh elements
  397. CALL CMISSMesh_NumberOfElementsSet(Mesh,TOTAL_NUMBER_OF_ELEMENTS,Err)
  398. !Set number of mesh components
  399. CALL CMISSMesh_NumberOfComponentsSet(Mesh,MESH_NUMBER_OF_COMPONENTS,Err)
  400. !Specify spatial mesh component
  401. CALL CMISSMeshElements_Initialise(MeshElementsSpace,Err)
  402. CALL CMISSMeshElements_Initialise(MeshElementsVelocity,Err)
  403. CALL CMISSMeshElements_Initialise(MeshElementsPressure,Err)
  404. MESH_COMPONENT_NUMBER_SPACE=1
  405. MESH_COMPONENT_NUMBER_VELOCITY=1
  406. MESH_COMPONENT_NUMBER_PRESSURE=1
  407. CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_SPACE,BasisSpace,MeshElementsSpace,Err)
  408. DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  409. CALL CMISSMeshElements_NodesSet(MeshElementsSpace,ELEMENT_NUMBER,CM%M(ELEMENT_NUMBER,1:NUMBER_OF_ELEMENT_NODES_SPACE),Err)
  410. ENDDO
  411. CALL CMISSMeshElements_CreateFinish(MeshElementsSpace,Err)
  412. !Specify velocity mesh component
  413. IF(BASIS_XI_INTERPOLATION_VELOCITY==BASIS_XI_INTERPOLATION_SPACE) THEN
  414. MeshElementsVelocity=MeshElementsSpace
  415. ELSE
  416. MESH_COMPONENT_NUMBER_VELOCITY=MESH_COMPONENT_NUMBER_SPACE+1
  417. CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_VELOCITY,BasisVelocity,MeshElementsVelocity,Err)
  418. DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  419. CALL CMISSMeshElements_NodesSet(MeshElementsVelocity,ELEMENT_NUMBER,CM%V(ELEMENT_NUMBER, &
  420. & 1:NUMBER_OF_ELEMENT_NODES_VELOCITY),Err)
  421. ENDDO
  422. CALL CMISSMeshElements_CreateFinish(MeshElementsVelocity,Err)
  423. ENDIF
  424. !Specify pressure mesh component
  425. IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_SPACE) THEN
  426. MeshElementsPressure=MeshElementsSpace
  427. MESH_COMPONENT_NUMBER_PRESSURE=MESH_COMPONENT_NUMBER_SPACE
  428. ELSE IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_VELOCITY) THEN
  429. MeshElementsPressure=MeshElementsVelocity
  430. MESH_COMPONENT_NUMBER_PRESSURE=MESH_COMPONENT_NUMBER_VELOCITY
  431. ELSE
  432. MESH_COMPONENT_NUMBER_PRESSURE=MESH_COMPONENT_NUMBER_VELOCITY+1
  433. CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_PRESSURE,BasisPressure,MeshElementsPressure,Err)
  434. DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  435. CALL CMISSMeshElements_NodesSet(MeshElementsPressure,ELEMENT_NUMBER,CM%P(ELEMENT_NUMBER, &
  436. & 1:NUMBER_OF_ELEMENT_NODES_PRESSURE),Err)
  437. ENDDO
  438. CALL CMISSMeshElements_CreateFinish(MeshElementsPressure,Err)
  439. ENDIF
  440. !Finish the creation of the mesh
  441. CALL CMISSMesh_CreateFinish(Mesh,Err)
  442. !
  443. !================================================================================================================================
  444. !
  445. !GEOMETRIC FIELD
  446. !Create a decomposition
  447. CALL CMISSDecomposition_Initialise(Decomposition,Err)
  448. CALL CMISSDecomposition_CreateStart(DecompositionUserNumber,Mesh,Decomposition,Err)
  449. !Set the decomposition to be a general decomposition with the specified number of domains
  450. CALL CMISSDecomposition_TypeSet(Decomposition,CMISS_DECOMPOSITION_CALCULATED_TYPE,Err)
  451. CALL CMISSDecomposition_NumberOfDomainsSet(Decomposition,NumberOfComputationalNodes,Err)
  452. !Finish the decomposition
  453. CALL CMISSDecomposition_CreateFinish(Decomposition,Err)
  454. !Start to create a default (geometric) field on the region
  455. CALL CMISSField_Initialise(GeometricField,Err)
  456. CALL CMISSField_CreateStart(GeometricFieldUserNumber,Region,GeometricField,Err)
  457. !Set the field type
  458. CALL CMISSField_TypeSet(GeometricField,CMISS_FIELD_GEOMETRIC_TYPE,Err)
  459. !Set the decomposition to use
  460. CALL CMISSField_MeshDecompositionSet(GeometricField,Decomposition,Err)
  461. !Set the scaling to use
  462. CALL CMISSField_ScalingTypeSet(GeometricField,CMISS_FIELD_NO_SCALING,Err)
  463. !Set the mesh component to be used by the field components.
  464. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  465. CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  466. & MESH_COMPONENT_NUMBER_SPACE,Err)
  467. ENDDO
  468. !Finish creating the field
  469. CALL CMISSField_CreateFinish(GeometricField,Err)
  470. !Update the geometric field parameters
  471. DO NODE_NUMBER=1,NUMBER_OF_NODES_SPACE
  472. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  473. VALUE=CM%N(NODE_NUMBER,COMPONENT_NUMBER)
  474. CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE_NUMBER,1,BoundaryNodeDomain,Err)
  475. IF(BoundaryNodeDomain==ComputationalNodeNumber) THEN
  476. CALL CMISSField_ParameterSetUpdateNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  477. & 1,CMISS_NO_GLOBAL_DERIV,NODE_NUMBER,COMPONENT_NUMBER,VALUE,Err)
  478. ENDIF
  479. ENDDO
  480. ENDDO
  481. CALL CMISSField_ParameterSetUpdateStart(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,Err)
  482. CALL CMISSField_ParameterSetUpdateFinish(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,Err)
  483. !
  484. !================================================================================================================================
  485. !
  486. !EQUATIONS SETS
  487. !Create the equations set for dynamic Stokes
  488. CALL CMISSEquationsSet_Initialise(EquationsSetStokes,Err)
  489. CALL CMISSField_Initialise(EquationsSetField,Err)
  490. CALL CMISSEquationsSet_CreateStart(EquationsSetUserNumberStokes,Region,GeometricField,CMISS_EQUATIONS_SET_FLUID_MECHANICS_CLASS, &
  491. & CMISS_EQUATIONS_SET_STOKES_EQUATION_TYPE,CMISS_EQUATIONS_SET_TRANSIENT_STOKES_SUBTYPE,EquationsSetFieldUserNumber, &
  492. & EquationsSetField, &
  493. & EquationsSetStokes,Err)
  494. !Set the equations set to be a dynamic Stokes problem
  495. !Finish creating the equations set
  496. CALL CMISSEquationsSet_CreateFinish(EquationsSetStokes,Err)
  497. !
  498. !================================================================================================================================
  499. !
  500. !DEPENDENT FIELDS
  501. !Create the equations set dependent field variables for dynamic Stokes
  502. CALL CMISSField_Initialise(DependentFieldStokes,Err)
  503. CALL CMISSEquationsSet_DependentCreateStart(EquationsSetStokes,DependentFieldUserNumberStokes, &
  504. & DependentFieldStokes,Err)
  505. !Set the mesh component to be used by the field components.
  506. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  507. CALL CMISSField_ComponentMeshComponentSet(DependentFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  508. & MESH_COMPONENT_NUMBER_VELOCITY,Err)
  509. CALL CMISSField_ComponentMeshComponentSet(DependentFieldStokes,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, &
  510. & MESH_COMPONENT_NUMBER_VELOCITY,Err)
  511. ENDDO
  512. COMPONENT_NUMBER=NUMBER_OF_DIMENSIONS+1
  513. CALL CMISSField_ComponentMeshComponentSet(DependentFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  514. & MESH_COMPONENT_NUMBER_PRESSURE,Err)
  515. CALL CMISSField_ComponentMeshComponentSet(DependentFieldStokes,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, &
  516. & MESH_COMPONENT_NUMBER_PRESSURE,Err)
  517. !Finish the equations set dependent field variables
  518. CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetStokes,Err)
  519. !Initialise dependent field
  520. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  521. CALL CMISSField_ComponentValuesInitialise(DependentFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  522. & COMPONENT_NUMBER,INITIAL_FIELD_STOKES(COMPONENT_NUMBER),Err)
  523. ENDDO
  524. !
  525. !================================================================================================================================
  526. !
  527. !MATERIALS FIELDS
  528. !Create the equations set materials field variables for dynamic Stokes
  529. CALL CMISSField_Initialise(MaterialsFieldStokes,Err)
  530. CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetStokes,MaterialsFieldUserNumberStokes, &
  531. & MaterialsFieldStokes,Err)
  532. !Finish the equations set materials field variables
  533. CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetStokes,Err)
  534. CALL CMISSField_ComponentValuesInitialise(MaterialsFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  535. & MaterialsFieldUserNumberStokesMu,MU_PARAM_STOKES,Err)
  536. CALL CMISSField_ComponentValuesInitialise(MaterialsFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  537. & MaterialsFieldUserNumberStokesRho,RHO_PARAM_STOKES,Err)
  538. !
  539. !================================================================================================================================
  540. !
  541. !EQUATIONS
  542. !Create the equations set equations
  543. CALL CMISSEquations_Initialise(EquationsStokes,Err)
  544. CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetStokes,EquationsStokes,Err)
  545. !Set the equations matrices sparsity type
  546. CALL CMISSEquations_SparsityTypeSet(EquationsStokes,CMISS_EQUATIONS_SPARSE_MATRICES,Err)
  547. !Set the equations lumping type
  548. CALL CMISSEquations_LumpingTypeSet(EquationsStokes,CMISS_EQUATIONS_UNLUMPED_MATRICES,Err)
  549. !Set the equations set output
  550. CALL CMISSEquations_OutputTypeSet(EquationsStokes,EQUATIONS_STOKES_OUTPUT,Err)
  551. !Finish the equations set equations
  552. CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetStokes,Err)
  553. !
  554. !================================================================================================================================
  555. !
  556. !PROBLEMS
  557. !Start the creation of a problem.
  558. CALL CMISSProblem_Initialise(Problem,Err)
  559. CALL CMISSControlLoop_Initialise(ControlLoop,Err)
  560. CALL CMISSProblem_CreateStart(ProblemUserNumber,Problem,Err)
  561. !Set the problem to be a dynamic Stokes problem
  562. CALL CMISSProblem_SpecificationSet(Problem,CMISS_PROBLEM_FLUID_MECHANICS_CLASS,CMISS_PROBLEM_STOKES_EQUATION_TYPE, &
  563. & CMISS_PROBLEM_TRANSIENT_STOKES_SUBTYPE,Err)
  564. !Finish the creation of a problem.
  565. CALL CMISSProblem_CreateFinish(Problem,Err)
  566. !Start the creation of the problem control loop
  567. CALL CMISSProblem_ControlLoopCreateStart(Problem,Err)
  568. !Get the control loop
  569. CALL CMISSProblem_ControlLoopGet(Problem,CMISS_CONTROL_LOOP_NODE,ControlLoop,Err)
  570. !Set the times
  571. CALL CMISSControlLoop_TimesSet(ControlLoop,DYNAMIC_SOLVER_STOKES_START_TIME,DYNAMIC_SOLVER_STOKES_STOP_TIME, &
  572. & DYNAMIC_SOLVER_STOKES_TIME_INCREMENT,Err)
  573. !Set the output timing
  574. CALL CMISSControlLoop_TimeOutputSet(ControlLoop,DYNAMIC_SOLVER_STOKES_OUTPUT_FREQUENCY,Err)
  575. !Finish creating the problem control loop
  576. CALL CMISSProblem_ControlLoopCreateFinish(Problem,Err)
  577. !
  578. !================================================================================================================================
  579. !
  580. !SOLVERS
  581. !Start the creation of the problem solvers
  582. CALL CMISSSolver_Initialise(DynamicSolverStokes,Err)
  583. CALL CMISSSolver_Initialise(LinearSolverStokes,Err)
  584. CALL CMISSProblem_SolversCreateStart(Problem,Err)
  585. !Get the dynamic dymamic solver
  586. CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,SolverStokesUserNumber,DynamicSolverStokes,Err)
  587. !Set the output type
  588. CALL CMISSSolver_OutputTypeSet(DynamicSolverStokes,DYNAMIC_SOLVER_STOKES_OUTPUT_TYPE,Err)
  589. !Set theta
  590. CALL CMISSSolver_DynamicThetaSet(DynamicSolverStokes,DYNAMIC_SOLVER_STOKES_THETA,Err)
  591. ! CALL CMISSSolverDynamicDynamicSet(DynamicSolverStokes,.TRUE.,Err)
  592. !Get the dynamic linear solver
  593. CALL CMISSSolver_DynamicLinearSolverGet(DynamicSolverStokes,LinearSolverStokes,Err)
  594. !Set the output type
  595. CALL CMISSSolver_OutputTypeSet(LinearSolverStokes,LINEAR_SOLVER_STOKES_OUTPUT_TYPE,Err)
  596. !Set the solver settings
  597. IF(LINEAR_SOLVER_STOKES_DIRECT_FLAG) THEN
  598. CALL CMISSSolver_LinearTypeSet(LinearSolverStokes,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err)
  599. CALL CMISSSolver_LibraryTypeSet(LinearSolverStokes,CMISS_SOLVER_MUMPS_LIBRARY,Err)
  600. ELSE
  601. CALL CMISSSolver_LinearTypeSet(LinearSolverStokes,CMISS_SOLVER_LINEAR_ITERATIVE_SOLVE_TYPE,Err)
  602. CALL CMISSSolver_LinearIterativeMaximumIterationsSet(LinearSolverStokes,MAXIMUM_ITERATIONS,Err)
  603. CALL CMISSSolver_LinearIterativeDivergenceToleranceSet(LinearSolverStokes,DIVERGENCE_TOLERANCE,Err)
  604. CALL CMISSSolver_LinearIterativeRelativeToleranceSet(LinearSolverStokes,RELATIVE_TOLERANCE,Err)
  605. CALL CMISSSolver_LinearIterativeAbsoluteToleranceSet(LinearSolverStokes,ABSOLUTE_TOLERANCE,Err)
  606. CALL CMISSSolver_LinearIterativeGMRESRestartSet(LinearSolverStokes,RESTART_VALUE,Err)
  607. ENDIF
  608. !Finish the creation of the problem solver
  609. CALL CMISSProblem_SolversCreateFinish(Problem,Err)
  610. !
  611. !================================================================================================================================
  612. !
  613. !SOLVER EQUATIONS
  614. !Start the creation of the problem solver equations
  615. CALL CMISSSolver_Initialise(DynamicSolverStokes,Err)
  616. CALL CMISSSolverEquations_Initialise(SolverEquationsStokes,Err)
  617. CALL CMISSProblem_SolverEquationsCreateStart(Problem,Err)
  618. !Get the dynamic solver equations
  619. CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,SolverStokesUserNumber,DynamicSolverStokes,Err)
  620. CALL CMISSSolver_SolverEquationsGet(DynamicSolverStokes,SolverEquationsStokes,Err)
  621. !Set the solver equations sparsity
  622. CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsStokes,CMISS_SOLVER_SPARSE_MATRICES,Err)
  623. !Add in the equations set
  624. CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsStokes,EquationsSetStokes,EquationsSetIndex,Err)
  625. !Finish the creation of the problem solver equations
  626. CALL CMISSProblem_SolverEquationsCreateFinish(Problem,Err)
  627. !
  628. !================================================================================================================================
  629. !
  630. !BOUNDARY CONDITIONS
  631. !Start the creation of the equations set boundary conditions for Stokes
  632. CALL CMISSBoundaryConditions_Initialise(BoundaryConditionsStokes,Err)
  633. CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsStokes,BoundaryConditionsStokes,Err)
  634. !Set fixed wall nodes
  635. IF(FIXED_WALL_NODES_STOKES_FLAG) THEN
  636. DO NODE_COUNTER=1,NUMBER_OF_FIXED_WALL_NODES_STOKES
  637. NODE_NUMBER=FIXED_WALL_NODES_STOKES(NODE_COUNTER)
  638. CONDITION=CMISS_BOUNDARY_CONDITION_FIXED_WALL
  639. CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE_NUMBER,1,BoundaryNodeDomain,Err)
  640. IF(BoundaryNodeDomain==ComputationalNodeNumber) THEN
  641. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  642. VALUE=0.0_CMISSDP
  643. CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsStokes,DependentFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,1, &
  644. & CMISS_NO_GLOBAL_DERIV, &
  645. & NODE_NUMBER,COMPONENT_NUMBER,CONDITION,VALUE,Err)
  646. ENDDO
  647. ENDIF
  648. ENDDO
  649. ENDIF
  650. !Set velocity boundary conditions
  651. IF(INLET_WALL_NODES_STOKES_FLAG) THEN
  652. DO NODE_COUNTER=1,NUMBER_OF_INLET_WALL_NODES_STOKES
  653. NODE_NUMBER=INLET_WALL_NODES_STOKES(NODE_COUNTER)
  654. CONDITION=CMISS_BOUNDARY_CONDITION_FIXED_INLET
  655. CALL CMISSDecomposition_NodeDomainGet(Decomposition,NODE_NUMBER,1,BoundaryNodeDomain,Err)
  656. IF(BoundaryNodeDomain==ComputationalNodeNumber) THEN
  657. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  658. VALUE=BOUNDARY_CONDITIONS_STOKES(COMPONENT_NUMBER)
  659. CALL CMISSBoundaryConditions_SetNode(BoundaryConditionsStokes,DependentFieldStokes,CMISS_FIELD_U_VARIABLE_TYPE,1, &
  660. & CMISS_NO_GLOBAL_DERIV, &
  661. & NODE_NUMBER,COMPONENT_NUMBER,CONDITION,VALUE,Err)
  662. ENDDO
  663. ENDIF
  664. ENDDO
  665. ENDIF
  666. !Finish the creation of the equations set boundary conditions
  667. CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquationsStokes,Err)
  668. !
  669. !================================================================================================================================
  670. !
  671. !RUN SOLVERS
  672. !Turn of PETSc error handling
  673. !CALL PETSC_ERRORHANDLING_SET_ON(ERR,ERROR,*999)
  674. !Solve the problem
  675. WRITE(*,'(A)') "Solving problem..."
  676. CALL CMISSProblem_Solve(Problem,Err)
  677. WRITE(*,'(A)') "Problem solved!"
  678. !
  679. !================================================================================================================================
  680. !
  681. !OUTPUT
  682. EXPORT_FIELD_IO=.TRUE.
  683. IF(EXPORT_FIELD_IO) THEN
  684. WRITE(*,'(A)') "Exporting fields..."
  685. CALL CMISSFields_Initialise(Fields,Err)
  686. CALL CMISSFields_Create(Region,Fields,Err)
  687. CALL CMISSFields_NodesExport(Fields,"DynamicStokes","FORTRAN",Err)
  688. CALL CMISSFields_ElementsExport(Fields,"DynamicStokes","FORTRAN",Err)
  689. CALL CMISSFields_Finalise(Fields,Err)
  690. WRITE(*,'(A)') "Field exported!"
  691. ENDIF
  692. !Finialise CMISS
  693. ! CALL CMISSFinalise(Err)
  694. WRITE(*,'(A)') "Program successfully completed."
  695. STOP
  696. END PROGRAM STOKESDYNAMICEXAMPLE