PageRenderTime 56ms CodeModel.GetById 17ms RepoModel.GetById 0ms app.codeStats 0ms

/MultiPhysics/CoupledDiffusionDiffusion/MonolithicSchemeTestAnalytic/src/MonolithicSchemeTestAnalyticExample.f90

http://github.com/xyan075/examples
FORTRAN Modern | 1129 lines | 516 code | 127 blank | 486 comment | 0 complexity | b0aa9e05b74e8f26ad63c016923ddb23 MD5 | raw file
  1. !> \file
  2. !> \authors Andrew Cookson
  3. !> \brief This is an example program to solve coupled multi-compartment diffusion equations in monolithic scheme 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 MultiPhysics/CoupledDiffusionDiffusion/MonolithicSchemeTest/src/MonolithicSchemeTestExample.f90
  42. !! Example program to solve coupled MonolithicSchemeTest equations using OpenCMISS calls.
  43. !! \par Latest Builds:
  44. !! \li <a href='http://autotest.bioeng.auckland.ac.nz/opencmiss-build/logs_x86_64-linux/MultiPhysics/CoupledDiffusionDiffusion/MonolithicSchemeTest/build-intel'>Linux Intel Build</a>
  45. !! \li <a href='http://autotest.bioeng.auckland.ac.nz/opencmiss-build/logs_x86_64-linux/MultiPhysics/CoupledDiffusionDiffusion/MonolithicSchemeTest/build-intel'>Linux GNU Build</a>
  46. !!
  47. !<
  48. ! !
  49. ! ! This example considers a volume coupled multi-compartment diffusion model, as a means of testing the monolithic assembly of
  50. ! ! such a system of equations, for use in more complicated problem in future.
  51. ! ! This example will initially couple together three diffusion equations, with the transfer between each equation being proportional to
  52. ! ! the concentration difference between each equation.
  53. !> Main program
  54. PROGRAM MONOLITHICSCHEMETESTFIELDMLEXAMPLE
  55. !
  56. !================================================================================================================================
  57. !
  58. !PROGRAM LIBRARIES
  59. USE OPENCMISS
  60. USE FIELDML_OUTPUT_ROUTINES
  61. USE FIELDML_UTIL_ROUTINES
  62. USE FIELDML_API
  63. USE FLUID_MECHANICS_IO_ROUTINES
  64. USE MPI
  65. #ifdef WIN32
  66. USE IFQWINCMISS
  67. #endif
  68. !
  69. !================================================================================================================================
  70. !
  71. !PROGRAM VARIABLES AND TYPES
  72. IMPLICIT NONE
  73. !Test program parameters
  74. REAL(CMISSDP), PARAMETER :: HEIGHT=1.0_CMISSDP
  75. REAL(CMISSDP), PARAMETER :: WIDTH=1.0_CMISSDP
  76. REAL(CMISSDP), PARAMETER :: LENGTH=3.0_CMISSDP
  77. INTEGER(CMISSIntg), PARAMETER :: CoordinateSystemUserNumber=1
  78. INTEGER(CMISSIntg), PARAMETER :: RegionUserNumber=2
  79. INTEGER(CMISSIntg), PARAMETER :: MeshUserNumber=3
  80. INTEGER(CMISSIntg), PARAMETER :: DecompositionUserNumber=4
  81. INTEGER(CMISSIntg), PARAMETER :: GeometricFieldUserNumber=5
  82. INTEGER(CMISSIntg), PARAMETER :: DependentFieldUserNumber=6
  83. INTEGER(CMISSIntg), PARAMETER :: GeneratedMeshUserNumber=9
  84. INTEGER(CMISSIntg) :: MaterialsFieldUserNumberDiffusion
  85. ! INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDiffusionTwo=10
  86. ! INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberDiffusionThree=11
  87. INTEGER(CMISSIntg) :: SourceFieldUserNumberDiffusion
  88. ! INTEGER(CMISSIntg), PARAMETER :: SourceFieldUserNumberDiffusionTwo=16
  89. ! INTEGER(CMISSIntg), PARAMETER :: SourceFieldUserNumberDiffusionThree=17
  90. INTEGER(CMISSIntg) :: EquationsSetUserNumberDiffusion
  91. INTEGER(CMISSIntg) :: AnalyticFieldUserNumberDiffusion
  92. ! INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumberDiffusionTwo=13
  93. ! INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumberDiffusionThree=14
  94. INTEGER(CMISSIntg), PARAMETER :: ProblemUserNumber=20
  95. !!!!!!!!!!!!!!!!!!!EQUATIONS FIELD PARAMETERS
  96. INTEGER(CMISSIntg), PARAMETER :: EquationsFieldDiffusionOne=33
  97. INTEGER(CMISSIntg), PARAMETER :: EquationsFieldDiffusionTwo=34
  98. INTEGER(CMISSIntg), PARAMETER :: EquationsFieldDiffusionThree=35
  99. INTEGER(CMISSIntg) :: EquationsSetFieldUserNumberDiffusion
  100. INTEGER(CMISSIntg) :: icompartment,Ncompartments,num_var
  101. !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
  102. INTEGER(CMISSIntg), PARAMETER :: DomainUserNumber=1
  103. INTEGER(CMISSIntg), PARAMETER :: NumberOfUserDomains=1
  104. INTEGER(CMISSIntg), PARAMETER :: SolverDiffusionUserNumber=1
  105. !Program types
  106. TYPE(EXPORT_CONTAINER):: CM
  107. !Program variables
  108. INTEGER(CMISSIntg) :: NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS,NUMBER_GLOBAL_Z_ELEMENTS
  109. INTEGER(CMISSIntg) :: NUMBER_OF_DOMAINS
  110. INTEGER(CMISSIntg) :: MPI_IERROR
  111. INTEGER(CMISSIntg) :: NUMBER_OF_DIMENSIONS
  112. INTEGER(CMISSIntg) :: NUMBER_OF_COMPARTMENTS
  113. INTEGER(CMISSIntg) :: BASIS_TYPE
  114. INTEGER(CMISSIntg) :: BASIS_NUMBER_GEOMETRY
  115. INTEGER(CMISSIntg) :: BASIS_NUMBER_CONC_ONE
  116. INTEGER(CMISSIntg) :: BASIS_NUMBER_CONC_TWO
  117. INTEGER(CMISSIntg) :: BASIS_NUMBER_CONC_THREE
  118. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_GEOMETRY
  119. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_CONC_ONE
  120. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_CONC_TWO
  121. INTEGER(CMISSIntg) :: BASIS_XI_GAUSS_CONC_THREE
  122. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_GEOMETRY
  123. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_CONC_ONE
  124. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_CONC_TWO
  125. INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_CONC_THREE
  126. INTEGER(CMISSIntg) :: MESH_NUMBER_OF_COMPONENTS,MESH_NUMBER_OF_ALL_COMPONENTS
  127. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_GEOMETRY
  128. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_CONC_ONE
  129. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_CONC_TWO
  130. INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_CONC_THREE
  131. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_GEOMETRY
  132. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_CONC_ONE
  133. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_CONC_TWO
  134. INTEGER(CMISSIntg) :: NUMBER_OF_NODES_CONC_THREE
  135. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_GEOMETRY
  136. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_CONC_ONE
  137. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_CONC_TWO
  138. INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_CONC_THREE
  139. INTEGER(CMISSIntg) :: TOTAL_NUMBER_OF_NODES,TOTAL_NUMBER_OF_ALL_NODES
  140. INTEGER(CMISSIntg) :: TOTAL_NUMBER_OF_ELEMENTS
  141. INTEGER(CMISSIntg) :: MAXIMUM_ITERATIONS
  142. INTEGER(CMISSIntg) :: RESTART_VALUE
  143. INTEGER(CMISSIntg) :: NUMBER_OF_FIXED_WALL_NODES_DIFFUSION_ONE
  144. INTEGER(CMISSIntg) :: NUMBER_OF_INLET_WALL_NODES_DIFFUSION_ONE
  145. INTEGER(CMISSIntg) :: NUMBER_OF_FIXED_WALL_NODES_DIFFUSION_TWO
  146. INTEGER(CMISSIntg) :: NUMBER_OF_INLET_WALL_NODES_DIFFUSION_TWO
  147. INTEGER(CMISSIntg) :: NUMBER_OF_FIXED_WALL_NODES_DIFFUSION_THREE
  148. INTEGER(CMISSIntg) :: NUMBER_OF_INLET_WALL_NODES_DIFFUSION_THREE
  149. INTEGER(CMISSIntg) :: EQUATIONS_DIFFUSION_OUTPUT
  150. INTEGER(CMISSIntg) :: EQUATIONS_DIFFUSION_TWO_OUTPUT
  151. INTEGER(CMISSIntg) :: EQUATIONS_DIFFUSION_THREE_OUTPUT
  152. INTEGER(CMISSIntg) :: COMPONENT_NUMBER
  153. INTEGER(CMISSIntg) :: NODE_NUMBER
  154. INTEGER(CMISSIntg) :: ELEMENT_NUMBER
  155. INTEGER(CMISSIntg) :: NODE_COUNTER
  156. INTEGER(CMISSIntg) :: CONDITION
  157. INTEGER(CMISSIntg) :: LINEAR_SOLVER_DIFFUSION_OUTPUT_FREQUENCY
  158. INTEGER(CMISSIntg) :: LINEAR_SOLVER_DIFFUSION_OUTPUT_TYPE
  159. REAL(CMISSDP) :: DOMAIN_X1, DOMAIN_X2, DOMAIN_Y1, DOMAIN_Y2, DOMAIN_Z1, DOMAIN_Z2
  160. REAL(CMISSDP) :: GEOMETRY_TOLERANCE
  161. INTEGER, ALLOCATABLE, DIMENSION(:):: FIXED_WALL_NODES_DIFFUSION_ONE
  162. INTEGER, ALLOCATABLE, DIMENSION(:):: INLET_WALL_NODES_DIFFUSION_ONE
  163. INTEGER, ALLOCATABLE, DIMENSION(:):: FIXED_WALL_NODES_DIFFUSION_TWO
  164. INTEGER, ALLOCATABLE, DIMENSION(:):: INLET_WALL_NODES_DIFFUSION_TWO
  165. INTEGER, ALLOCATABLE, DIMENSION(:):: FIXED_WALL_NODES_DIFFUSION_THREE
  166. INTEGER, ALLOCATABLE, DIMENSION(:):: INLET_WALL_NODES_DIFFUSION_THREE
  167. REAL(CMISSDP) :: INITIAL_FIELD_DIFFUSION_ONE
  168. REAL(CMISSDP) :: INITIAL_FIELD_DIFFUSION_TWO
  169. REAL(CMISSDP) :: INITIAL_FIELD_DIFFUSION_THREE
  170. REAL(CMISSDP) :: BOUNDARY_CONDITIONS_DIFFUSION_ONE
  171. REAL(CMISSDP) :: BOUNDARY_CONDITIONS_DIFFUSION_TWO
  172. REAL(CMISSDP) :: BOUNDARY_CONDITIONS_DIFFUSION_THREE
  173. REAL(CMISSDP) :: DIVERGENCE_TOLERANCE
  174. REAL(CMISSDP) :: RELATIVE_TOLERANCE
  175. REAL(CMISSDP) :: ABSOLUTE_TOLERANCE
  176. REAL(CMISSDP) :: LINESEARCH_ALPHA
  177. REAL(CMISSDP) :: VALUE
  178. REAL(CMISSDP) :: LINEAR_SOLVER_DIFFUSION_START_TIME
  179. REAL(CMISSDP) :: LINEAR_SOLVER_DIFFUSION_STOP_TIME
  180. REAL(CMISSDP) :: LINEAR_SOLVER_DIFFUSION_TIME_INCREMENT
  181. LOGICAL :: EXPORT_FIELD_IO
  182. LOGICAL :: LINEAR_SOLVER_DIFFUSION_DIRECT_FLAG
  183. LOGICAL :: INLET_WALL_NODES_DIFFUSION_ONE_FLAG
  184. LOGICAL :: INLET_WALL_NODES_DIFFUSION_TWO_FLAG
  185. LOGICAL :: INLET_WALL_NODES_DIFFUSION_THREE_FLAG
  186. !CMISS variables
  187. CHARACTER(C_CHAR), PARAMETER :: NUL=C_NULL_CHAR
  188. !Regions
  189. TYPE(CMISSRegionType) :: Region
  190. TYPE(CMISSRegionType) :: WorldRegion
  191. !Coordinate systems
  192. TYPE(CMISSCoordinateSystemType) :: CoordinateSystem
  193. TYPE(CMISSCoordinateSystemType) :: WorldCoordinateSystem
  194. !Basis
  195. TYPE(CMISSBasisType) :: BasisGeometry
  196. TYPE(CMISSBasisType) :: BasisConcOne
  197. TYPE(CMISSBasisType) :: BasisConcTwo
  198. TYPE(CMISSBasisType) :: BasisConcThree
  199. !TYPE(CMISSBasisType), ALLOCATABLE, DIMENSION(:) :: BasicConc
  200. !Nodes
  201. TYPE(CMISSNodesType) :: Nodes
  202. !Elements
  203. TYPE(CMISSMeshElementsType) :: MeshElementsGeometry
  204. TYPE(CMISSMeshElementsType) :: MeshElementsConcOne
  205. TYPE(CMISSMeshElementsType) :: MeshElementsConcTwo
  206. TYPE(CMISSMeshElementsType) :: MeshElementsConcThree
  207. TYPE(CMISSMeshElementsType), ALLOCATABLE, DIMENSION(:) :: MeshElementsConc
  208. !Meshes
  209. TYPE(CMISSMeshType) :: Mesh
  210. !Decompositions
  211. TYPE(CMISSDecompositionType) :: Decomposition
  212. !Fields
  213. TYPE(CMISSFieldsType) :: Fields
  214. TYPE(CMISSGeneratedMeshType) :: GeneratedMesh
  215. !Field types
  216. TYPE(CMISSFieldType) :: GeometricField
  217. TYPE(CMISSFieldType) :: DependentField
  218. TYPE(CMISSFieldType), ALLOCATABLE, DIMENSION(:) :: MaterialsFieldDiffusion
  219. TYPE(CMISSFieldType), ALLOCATABLE, DIMENSION(:) :: SourceFieldDiffusion
  220. TYPE(CMISSFieldType), ALLOCATABLE, DIMENSION(:) :: AnalyticFieldDiffusion
  221. TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsDiffusion
  222. TYPE(CMISSEquationsSetType), ALLOCATABLE, DIMENSION(:) :: EquationsSetDiffusion
  223. TYPE(CMISSEquationsType), ALLOCATABLE, DIMENSION(:) :: EquationsDiffusion
  224. TYPE(CMISSFieldType), ALLOCATABLE, DIMENSION(:) :: EquationsSetFieldDiffusion
  225. !Problems
  226. TYPE(CMISSProblemType) :: Problem
  227. !Control loops
  228. TYPE(CMISSControlLoopType) :: ControlLoop
  229. !Solvers
  230. TYPE(CMISSSolverType) :: SolverDiffusion
  231. TYPE(CMISSSolverType) :: LinearSolverDiffusion
  232. !Solver equations
  233. TYPE(CMISSSolverEquationsType) :: SolverEquationsDiffusion
  234. #ifdef WIN32
  235. !Quickwin type
  236. LOGICAL :: QUICKWIN_STATUS=.FALSE.
  237. TYPE(WINDOWCONFIG) :: QUICKWIN_WINDOW_CONFIG
  238. #endif
  239. !Generic CMISS variables
  240. INTEGER(CMISSIntg) :: NumberOfComputationalNodes,ComputationalNodeNumber
  241. INTEGER(CMISSIntg) :: EquationsSetIndex
  242. INTEGER(CMISSIntg) :: EquationsSetIndexTwo
  243. INTEGER(CMISSIntg) :: EquationsSetIndexThree
  244. INTEGER(CMISSIntg) :: Err
  245. !Array containing the field variable types that will be used (for ease of incorporating inside a loop)
  246. INTEGER(CMISSIntg), ALLOCATABLE, DIMENSION(:) :: VariableTypes
  247. REAL(CMISSDP), ALLOCATABLE, DIMENSION(:,:) :: CouplingCoeffs
  248. INTEGER(CMISSIntg) :: DIAG_LEVEL_LIST(5)
  249. ! CHARACTER(LEN=255) :: DIAG_ROUTINE_LIST(8) !,TIMING_ROUTINE_LIST(1)
  250. CHARACTER(LEN=255) :: DIAG_ROUTINE_LIST(1) !,TIMING_ROUTINE_LIST(1)
  251. INTEGER(CMISSIntg) :: TotalNumberOfSolidNodes
  252. ! INTEGER(CMISSIntg) :: NumberOfSolidMeshComponents
  253. !
  254. !--------------------------------------------------------------------------------------------------------------------------------
  255. !
  256. !FieldML variables
  257. CHARACTER(KIND=C_CHAR,LEN=*), PARAMETER :: outputDirectory = ""
  258. CHARACTER(KIND=C_CHAR,LEN=*), PARAMETER :: outputFilename = "MonolithicMultiCompDiffusion.xml"
  259. CHARACTER(KIND=C_CHAR,LEN=*), PARAMETER :: basename = "monolithic_multicomp_diffusion"
  260. TYPE(FieldmlInfoType) :: fieldmlInfo
  261. #ifdef WIN32
  262. !Initialise QuickWin
  263. QUICKWIN_WINDOW_CONFIG%TITLE="General Output" !Window title
  264. QUICKWIN_WINDOW_CONFIG%NUMTEXTROWS=-1 !Max possible number of rows
  265. QUICKWIN_WINDOW_CONFIG%MODE=QWIN$SCROLLDOWN
  266. !Set the window parameters
  267. QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
  268. !If attempt fails set with system estimated values
  269. IF(.NOT.QUICKWIN_STATUS) QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
  270. #endif
  271. !
  272. !================================================================================================================================
  273. !
  274. !INITIALISE OPENCMISS
  275. CALL CMISSInitialise(WorldCoordinateSystem,WorldRegion,Err)
  276. !CALL CMISSErrorHandlingModeSet(CMISS_ERRORS_TRAP_ERROR,Err)
  277. !
  278. !================================================================================================================================
  279. !
  280. !PROBLEM CONTROL PANEL
  281. !Import cmHeart mesh information
  282. CALL FLUID_MECHANICS_IO_READ_CMHEART(CM,Err)
  283. BASIS_NUMBER_GEOMETRY=CM%ID_M
  284. BASIS_NUMBER_CONC_ONE=CM%ID_V !USE THE V COMPONENT OF CMHEART INPUT FOR CONCENTRATION ONE
  285. BASIS_NUMBER_CONC_TWO=CM%ID_P !USE THE p COMPONENT OF CMHEART INPUT FOR CONCENTRATION TWO
  286. NUMBER_OF_DIMENSIONS=CM%D
  287. BASIS_TYPE=CM%IT_T
  288. BASIS_XI_INTERPOLATION_GEOMETRY=CM%IT_M
  289. ! BASIS_XI_INTERPOLATION_GEOMETRY=2
  290. BASIS_XI_INTERPOLATION_CONC_ONE=CM%IT_V
  291. BASIS_XI_INTERPOLATION_CONC_TWO=CM%IT_P
  292. NUMBER_OF_NODES_GEOMETRY=CM%N_M
  293. NUMBER_OF_NODES_CONC_ONE=CM%N_V
  294. NUMBER_OF_NODES_CONC_TWO=CM%N_P
  295. TOTAL_NUMBER_OF_NODES=CM%N_T
  296. TOTAL_NUMBER_OF_ELEMENTS=CM%E_T
  297. NUMBER_OF_ELEMENT_NODES_GEOMETRY=CM%EN_M
  298. NUMBER_OF_ELEMENT_NODES_CONC_ONE=CM%EN_V
  299. NUMBER_OF_ELEMENT_NODES_CONC_TWO=CM%EN_P
  300. ! !Set domain dimensions
  301. ! DOMAIN_X1 = 0.0_CMISSDP
  302. ! DOMAIN_X2 = 1.0_CMISSDP
  303. ! DOMAIN_Y1 = 0.0_CMISSDP
  304. ! DOMAIN_Y2 = 1.0_CMISSDP
  305. ! DOMAIN_Z1 = 0.0_CMISSDP
  306. ! DOMAIN_Z2 = 1.0_CMISSDP
  307. ! !Set geometric tolerance
  308. ! GEOMETRY_TOLERANCE = 1.0E-12_CMISSDP
  309. ! !Set initial values
  310. ! INITIAL_FIELD_DIFFUSION_ONE=1.0_CMISSDP
  311. ! INITIAL_FIELD_DIFFUSION_TWO=1.0_CMISSDP
  312. ! INITIAL_FIELD_DIFFUSION_THREE=1.0_CMISSDP
  313. ! !Set initial boundary conditions
  314. ! INLET_WALL_NODES_DIFFUSION_ONE_FLAG=.TRUE.
  315. ! IF(INLET_WALL_NODES_DIFFUSION_ONE_FLAG) THEN
  316. ! NUMBER_OF_INLET_WALL_NODES_DIFFUSION_ONE=36
  317. ! ALLOCATE(INLET_WALL_NODES_DIFFUSION_ONE(NUMBER_OF_INLET_WALL_NODES_DIFFUSION_ONE))
  318. ! INLET_WALL_NODES_DIFFUSION_ONE=(/191,155,119,83,23,21,192,156,120,84,24,22,&
  319. ! & 198,162,126,90,36,35,204,168,132,96,48,47,210,174,138,102,60,59,216,180,144,108,72,71/)
  320. ! !Set initial boundary conditions
  321. ! BOUNDARY_CONDITIONS_DIFFUSION_ONE=1.0_CMISSDP
  322. ! ENDIF !Set material parameters
  323. !
  324. ! INLET_WALL_NODES_DIFFUSION_TWO_FLAG=.TRUE.
  325. ! IF(INLET_WALL_NODES_DIFFUSION_TWO_FLAG) THEN
  326. ! NUMBER_OF_INLET_WALL_NODES_DIFFUSION_TWO=36
  327. ! ALLOCATE(INLET_WALL_NODES_DIFFUSION_TWO(NUMBER_OF_INLET_WALL_NODES_DIFFUSION_TWO))
  328. ! INLET_WALL_NODES_DIFFUSION_TWO=(/191,155,119,83,23,21,192,156,120,84,24,22,&
  329. ! & 198,162,126,90,36,35,204,168,132,96,48,47,210,174,138,102,60,59,216,180,144,108,72,71/)
  330. ! !Set initial boundary conditions
  331. ! BOUNDARY_CONDITIONS_DIFFUSION_TWO=1.0_CMISSDP
  332. ! ENDIF !Set material parameters
  333. !
  334. ! INLET_WALL_NODES_DIFFUSION_THREE_FLAG=.TRUE.
  335. ! IF(INLET_WALL_NODES_DIFFUSION_THREE_FLAG) THEN
  336. ! NUMBER_OF_INLET_WALL_NODES_DIFFUSION_THREE=36
  337. ! ALLOCATE(INLET_WALL_NODES_DIFFUSION_THREE(NUMBER_OF_INLET_WALL_NODES_DIFFUSION_THREE))
  338. ! INLET_WALL_NODES_DIFFUSION_THREE=(/191,155,119,83,23,21,192,156,120,84,24,22,&
  339. ! & 198,162,126,90,36,35,204,168,132,96,48,47,210,174,138,102,60,59,216,180,144,108,72,71/)
  340. ! !Set initial boundary conditions
  341. ! BOUNDARY_CONDITIONS_DIFFUSION_THREE=1.0_CMISSDP
  342. ! ENDIF !Set material parameters
  343. !Get the computational nodes information
  344. CALL CMISSComputationalNumberOfNodesGet(NumberOfComputationalNodes,Err)
  345. CALL CMISSComputationalNodeNumberGet(ComputationalNodeNumber,Err)
  346. NUMBER_GLOBAL_X_ELEMENTS=40
  347. NUMBER_GLOBAL_Y_ELEMENTS=40
  348. NUMBER_GLOBAL_Z_ELEMENTS=0
  349. IF(NUMBER_GLOBAL_Z_ELEMENTS==0)THEN
  350. NUMBER_OF_DIMENSIONS=2
  351. ELSE
  352. NUMBER_OF_DIMENSIONS=3
  353. ENDIF
  354. NUMBER_OF_DOMAINS=NumberOfComputationalNodes
  355. !Set all diganostic levels on for testing
  356. CALL MPI_BCAST(NUMBER_GLOBAL_X_ELEMENTS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR)
  357. CALL MPI_BCAST(NUMBER_GLOBAL_Y_ELEMENTS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR)
  358. CALL MPI_BCAST(NUMBER_GLOBAL_Z_ELEMENTS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR)
  359. CALL MPI_BCAST(NUMBER_OF_DOMAINS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR)
  360. !Set material parameters
  361. !Set number of Gauss points (Mind that also material field may be interpolated)
  362. BASIS_XI_GAUSS_GEOMETRY=3 !4
  363. BASIS_XI_GAUSS_CONC_ONE=3 !4
  364. BASIS_XI_GAUSS_CONC_TWO=3 !4
  365. !Set output parameter
  366. !(NoOutput/ProgressOutput/TimingOutput/SolverOutput/SolverMatrixOutput)
  367. LINEAR_SOLVER_DIFFUSION_OUTPUT_TYPE=CMISS_SOLVER_MATRIX_OUTPUT
  368. !LINEAR_SOLVER_DIFFUSION_OUTPUT_TYPE=CMISS_SOLVER_NO_OUTPUT
  369. !(NoOutput/TimingOutput/MatrixOutput/ElementOutput)
  370. EQUATIONS_DIFFUSION_OUTPUT=CMISS_EQUATIONS_MATRIX_OUTPUT
  371. !EQUATIONS_DIFFUSION_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT
  372. EQUATIONS_DIFFUSION_TWO_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT
  373. EQUATIONS_DIFFUSION_THREE_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT
  374. !Set time parameter
  375. LINEAR_SOLVER_DIFFUSION_START_TIME=0.0_CMISSDP
  376. LINEAR_SOLVER_DIFFUSION_STOP_TIME=0.3000001_CMISSDP
  377. LINEAR_SOLVER_DIFFUSION_TIME_INCREMENT=0.005_CMISSDP
  378. !Set result output parameter
  379. LINEAR_SOLVER_DIFFUSION_OUTPUT_FREQUENCY=1
  380. !Set solver parameters
  381. LINEAR_SOLVER_DIFFUSION_DIRECT_FLAG=.FALSE.
  382. RELATIVE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-05_CMISSDP
  383. ABSOLUTE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-10_CMISSDP
  384. DIVERGENCE_TOLERANCE=1.0E5_CMISSDP !default: 1.0E5
  385. MAXIMUM_ITERATIONS=10000_CMISSIntg !default: 100000
  386. RESTART_VALUE=30_CMISSIntg !default: 30
  387. LINESEARCH_ALPHA=1.0_CMISSDP
  388. icompartment =1_CMISSIntg
  389. Ncompartments=2_CMISSIntg
  390. !
  391. !================================================================================================================================
  392. !
  393. !Set diagnostics
  394. !DIAG_LEVEL_LIST(1)=1
  395. ! DIAG_LEVEL_LIST(2)=2
  396. ! DIAG_LEVEL_LIST(3)=3
  397. ! DIAG_LEVEL_LIST(4)=4
  398. ! DIAG_LEVEL_LIST(5)=5
  399. ! DIAG_ROUTINE_LIST(1)="DIFFUSION_EQUATION_FINITE_ELEMENT_CALCULATE"
  400. ! DIAG_ROUTINE_LIST(2)="PROBLEM_SOLVER_EQUATIONS_SOLVE"
  401. !DIAG_ROUTINE_LIST(1)="SOLVER_DYNAMIC_CREATE_FINISH"
  402. !CMISS_ALL_DIAG_TYPE/CMISS_IN_DIAG_TYPE/CMISS_FROM_DIAG_TYPE
  403. !CALL CMISSDiagnosticsSetOn(CMISS_IN_DIAG_TYPE,DIAG_LEVEL_LIST,"Diagnostics",DIAG_ROUTINE_LIST,Err)
  404. !CMISS_ALL_TIMING_TYPE/CMISS_IN_TIMING_TYPE/CMISS_FROM_TIMING_TYPE
  405. !TIMING_ROUTINE_LIST(1)="PROBLEM_FINITE_ELEMENT_CALCULATE"
  406. !CALL TIMING_SET_ON(IN_TIMING_TYPE,.TRUE.,"",TIMING_ROUTINE_LIST,ERR,ERROR,*999)
  407. !ALLOCATE THE ARRAYS
  408. ALLOCATE (EquationsSetDiffusion(Ncompartments))
  409. ALLOCATE (EquationsSetFieldDiffusion(Ncompartments))
  410. ALLOCATE (MaterialsFieldDiffusion(Ncompartments))
  411. ALLOCATE (SourceFieldDiffusion(Ncompartments))
  412. ALLOCATE (EquationsDiffusion(Ncompartments))
  413. ALLOCATE (AnalyticFieldDiffusion(Ncompartments))
  414. !
  415. !================================================================================================================================
  416. !
  417. !COORDINATE SYSTEM
  418. !Start the creation of a new RC coordinate system
  419. CALL CMISSCoordinateSystem_Initialise(CoordinateSystem,Err)
  420. CALL CMISSCoordinateSystem_CreateStart(CoordinateSystemUserNumber,CoordinateSystem,Err)
  421. !Set the coordinate system dimension
  422. CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,NUMBER_OF_DIMENSIONS,Err)
  423. !Finish the creation of the coordinate system
  424. CALL CMISSCoordinateSystem_CreateFinish(CoordinateSystem,Err)
  425. !
  426. !================================================================================================================================
  427. !
  428. !REGION
  429. !For a volume-coupled problem, both concentrations are based in the same region
  430. !Start the creation of a new region
  431. CALL CMISSRegion_Initialise(Region,Err)
  432. CALL CMISSRegion_CreateStart(RegionUserNumber,WorldRegion,Region,Err)
  433. !Set the regions coordinate system as defined above
  434. CALL CMISSRegion_CoordinateSystemSet(Region,CoordinateSystem,Err)
  435. !Finish the creation of the region
  436. CALL CMISSRegion_CreateFinish(Region,Err)
  437. !
  438. !================================================================================================================================
  439. !
  440. !BASES
  441. !Start the creation of new bases: Geometry
  442. MESH_NUMBER_OF_COMPONENTS=1
  443. CALL CMISSBasis_Initialise(BasisGeometry,Err)
  444. CALL CMISSBasis_CreateStart(BASIS_NUMBER_GEOMETRY,BasisGeometry,Err)
  445. !Set the basis type (Lagrange/Simplex)
  446. CALL CMISSBasis_TypeSet(BasisGeometry,BASIS_TYPE,Err)
  447. !Set the basis xi number
  448. CALL CMISSBasis_NumberOfXiSet(BasisGeometry,NUMBER_OF_DIMENSIONS,Err)
  449. !Set the basis xi interpolation and number of Gauss points
  450. IF(NUMBER_OF_DIMENSIONS==2) THEN
  451. CALL CMISSBasis_InterpolationXiSet(BasisGeometry,(/BASIS_XI_INTERPOLATION_GEOMETRY,BASIS_XI_INTERPOLATION_GEOMETRY/),Err)
  452. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisGeometry,(/BASIS_XI_GAUSS_GEOMETRY,BASIS_XI_GAUSS_GEOMETRY/),Err)
  453. ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
  454. CALL CMISSBasis_InterpolationXiSet(BasisGeometry,(/BASIS_XI_INTERPOLATION_GEOMETRY,BASIS_XI_INTERPOLATION_GEOMETRY, &
  455. & BASIS_XI_INTERPOLATION_GEOMETRY/),Err)
  456. CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisGeometry,(/BASIS_XI_GAUSS_GEOMETRY,BASIS_XI_GAUSS_GEOMETRY, &
  457. & BASIS_XI_GAUSS_GEOMETRY/),Err)
  458. ENDIF
  459. !Finish the creation of the basis
  460. CALL CMISSBasis_CreateFinish(BasisGeometry,Err)
  461. !
  462. !Start the creation of another basis: Concentration_One
  463. ! IF(BASIS_XI_INTERPOLATION_CONC_ONE==BASIS_XI_INTERPOLATION_GEOMETRY) THEN
  464. BasisConcOne=BasisGeometry
  465. ! ELSE
  466. ! MESH_NUMBER_OF_COMPONENTS=MESH_NUMBER_OF_COMPONENTS+1
  467. ! !Initialise a new velocity basis
  468. ! CALL CMISSBasis_Initialise(BasisConcOne,Err)
  469. ! !Start the creation of a basis
  470. ! CALL CMISSBasis_CreateStart(BASIS_NUMBER_CONC_ONE,BasisConcOne,Err)
  471. ! !Set the basis type (Lagrange/Simplex)
  472. ! CALL CMISSBasis_TypeSet(BasisConcOne,BASIS_TYPE,Err)
  473. ! !Set the basis xi number
  474. ! CALL CMISSBasis_NumberOfXiSet(BasisConcOne,NUMBER_OF_DIMENSIONS,Err)
  475. ! !Set the basis xi interpolation and number of Gauss points
  476. ! IF(NUMBER_OF_DIMENSIONS==2) THEN
  477. ! CALL CMISSBasis_InterpolationXiSet(BasisConcOne,(/BASIS_XI_INTERPOLATION_CONC_ONE,BASIS_XI_INTERPOLATION_CONC_ONE/),Err)
  478. ! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisConcOne,(/BASIS_XI_GAUSS_CONC_ONE,BASIS_XI_GAUSS_CONC_ONE/),Err)
  479. ! ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
  480. ! CALL CMISSBasis_InterpolationXiSet(BasisConcOne,(/BASIS_XI_INTERPOLATION_CONC_ONE,BASIS_XI_INTERPOLATION_CONC_ONE, &
  481. ! & BASIS_XI_INTERPOLATION_CONC_ONE/),Err)
  482. ! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisConcOne,(/BASIS_XI_GAUSS_CONC_ONE,BASIS_XI_GAUSS_CONC_ONE, &
  483. ! & BASIS_XI_GAUSS_CONC_ONE/),Err)
  484. ! ENDIF
  485. ! !Finish the creation of the basis
  486. ! CALL CMISSBasis_CreateFinish(BasisConcOne,Err)
  487. ! ENDIF
  488. !
  489. !Start the creation of another basis: Concentration_Two
  490. ! IF(BASIS_XI_INTERPOLATION_CONC_TWO==BASIS_XI_INTERPOLATION_GEOMETRY) THEN
  491. BasisConcTwo=BasisGeometry
  492. ! ELSE IF(BASIS_XI_INTERPOLATION_CONC_TWO==BASIS_XI_INTERPOLATION_CONC_ONE) THEN
  493. ! BasisConcTwo=BasisConcOne
  494. ! ELSE
  495. ! MESH_NUMBER_OF_COMPONENTS=MESH_NUMBER_OF_COMPONENTS+1
  496. ! !Initialise a new concentration basis
  497. ! CALL CMISSBasis_Initialise(BasisConcTwo,Err)
  498. ! !Start the creation of a basis
  499. ! CALL CMISSBasis_CreateStart(BASIS_NUMBER_CONC_TWO,BasisConcTwo,Err)
  500. ! !Set the basis type (Lagrange/Simplex)
  501. ! CALL CMISSBasis_TypeSet(BasisConcTwo,BASIS_TYPE,Err)
  502. ! !Set the basis xi number
  503. ! CALL CMISSBasis_NumberOfXiSet(BasisConcTwo,NUMBER_OF_DIMENSIONS,Err)
  504. ! !Set the basis xi interpolation and number of Gauss points
  505. ! IF(NUMBER_OF_DIMENSIONS==2) THEN
  506. ! CALL CMISSBasis_InterpolationXiSet(BasisConcTwo,(/BASIS_XI_INTERPOLATION_CONC_TWO,BASIS_XI_INTERPOLATION_CONC_TWO/),Err)
  507. ! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisConcTwo,(/BASIS_XI_GAUSS_CONC_TWO,BASIS_XI_GAUSS_CONC_TWO/),Err)
  508. ! ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
  509. ! CALL CMISSBasis_InterpolationXiSet(BasisConcTwo,(/BASIS_XI_INTERPOLATION_CONC_TWO,BASIS_XI_INTERPOLATION_CONC_TWO, &
  510. ! & BASIS_XI_INTERPOLATION_CONC_TWO/),Err)
  511. ! CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisConcTwo,(/BASIS_XI_GAUSS_CONC_TWO,BASIS_XI_GAUSS_CONC_TWO, &
  512. ! & BASIS_XI_GAUSS_CONC_TWO/),Err)
  513. ! ENDIF
  514. ! !Finish the creation of the basis
  515. ! CALL CMISSBasis_CreateFinish(BasisConcTwo,Err)
  516. ! ENDIF
  517. BasisConcThree=BasisGeometry
  518. !
  519. !================================================================================================================================
  520. !
  521. ! !MESH
  522. ! !All types of physics utilize the same "mesh", but may be represented on individual mesh components.
  523. ! TotalNumberOfSolidNodes = NUMBER_OF_NODES_GEOMETRY
  524. ! TOTAL_NUMBER_OF_ALL_NODES = TOTAL_NUMBER_OF_NODES + TotalNumberOfSolidNodes
  525. ! MESH_NUMBER_OF_ALL_COMPONENTS = MESH_NUMBER_OF_COMPONENTS
  526. !
  527. ! !Start the creation of mesh nodes
  528. ! CALL CMISSNodes_Initialise(Nodes,Err)
  529. ! CALL CMISSNodes_CreateStart(Region,TOTAL_NUMBER_OF_ALL_NODES,Nodes,Err)
  530. ! CALL CMISSNodes_CreateFinish(Nodes,Err)
  531. ! !Start the creation of the mesh
  532. ! CALL CMISSMesh_Initialise(Mesh,Err)
  533. ! CALL CMISSMesh_CreateStart(MeshUserNumber,Region,NUMBER_OF_DIMENSIONS,Mesh,Err)
  534. ! !Set number of mesh elements
  535. ! CALL CMISSMesh_NumberOfElementsSet(Mesh,TOTAL_NUMBER_OF_ELEMENTS,Err)
  536. ! !Set number of mesh components
  537. ! CALL CMISSMesh_NumberOfComponentsSet(Mesh,MESH_NUMBER_OF_ALL_COMPONENTS,Err)
  538. ! !
  539. ! CALL CMISSMeshElements_Initialise(MeshElementsGeometry,Err)
  540. ! CALL CMISSMeshElements_Initialise(MeshElementsConcOne,Err)
  541. ! CALL CMISSMeshElements_Initialise(MeshElementsConcTwo,Err)
  542. ! CALL CMISSMeshElements_Initialise(MeshElementsConcThree,Err)
  543. MESH_COMPONENT_NUMBER_GEOMETRY=1
  544. MESH_COMPONENT_NUMBER_CONC_ONE=1
  545. MESH_COMPONENT_NUMBER_CONC_TWO=1
  546. MESH_COMPONENT_NUMBER_CONC_THREE=1
  547. ! !Specify spatial mesh component
  548. ! CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_GEOMETRY,BasisGeometry,MeshElementsGeometry,Err)
  549. ! DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  550. ! CALL CMISSMeshElements_NodesSet(MeshElementsGeometry,ELEMENT_NUMBER,CM%M(ELEMENT_NUMBER,1:NUMBER_OF_ELEMENT_NODES_GEOMETRY),Err)
  551. ! ENDDO
  552. ! CALL CMISSMeshElements_CreateFinish(MeshElementsGeometry,Err)
  553. ! !Specify concentration one mesh component
  554. ! ! ! IF(BASIS_XI_INTERPOLATION_CONC_ONE==BASIS_XI_INTERPOLATION_GEOMETRY) THEN
  555. ! MeshElementsConcOne=MeshElementsGeometry
  556. ! ! ! ELSE
  557. ! ! ! MESH_COMPONENT_NUMBER_CONC_ONE=MESH_COMPONENT_NUMBER_GEOMETRY+1
  558. ! ! CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_CONC_ONE,BasisConcOne,MeshElementsConcOne,Err)
  559. ! ! ! DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  560. ! ! ! CALL CMISSMeshElements_NodesSet(MeshElementsConcOne,ELEMENT_NUMBER,CM%V(ELEMENT_NUMBER, &
  561. ! ! ! & 1:NUMBER_OF_ELEMENT_NODES_CONC_ONE),Err)
  562. ! ! ! ENDDO
  563. ! ! CALL CMISSMeshElements_CreateFinish(MeshElementsConcOne,Err)
  564. ! ! ! ENDIF
  565. ! ! !Specify concentration two mesh component
  566. ! ! ! IF(BASIS_XI_INTERPOLATION_CONC_TWO==BASIS_XI_INTERPOLATION_GEOMETRY) THEN
  567. ! MeshElementsConcTwo=MeshElementsGeometry
  568. ! ! ! MESH_COMPONENT_NUMBER_CONC_TWO=MESH_COMPONENT_NUMBER_GEOMETRY
  569. ! ! ! ELSE IF(BASIS_XI_INTERPOLATION_CONC_TWO==BASIS_XI_INTERPOLATION_CONC_ONE) THEN
  570. ! ! ! MeshElementsConcTwo=MeshElementsConcOne
  571. ! ! ! MESH_COMPONENT_NUMBER_CONC_TWO=MESH_COMPONENT_NUMBER_CONC_ONE
  572. ! ! ! ELSE
  573. ! ! ! MESH_COMPONENT_NUMBER_CONC_TWO=MESH_COMPONENT_NUMBER_CONC_ONE+1
  574. ! ! CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_CONC_TWO,BasisConcTwo,MeshElementsConcTwo,Err)
  575. ! ! ! DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  576. ! ! ! CALL CMISSMeshElements_NodesSet(MeshElementsConcTwo,ELEMENT_NUMBER,CM%P(ELEMENT_NUMBER, &
  577. ! ! ! & 1:NUMBER_OF_ELEMENT_NODES_CONC_TWO),Err)
  578. ! ! ! ENDDO
  579. ! ! CALL CMISSMeshElements_CreateFinish(MeshElementsConcTwo,Err)
  580. ! ! ! ENDIF
  581. ! ! !Specify concentration three mesh component
  582. ! ! ! IF(BASIS_XI_INTERPOLATION_CONC_THREE==BASIS_XI_INTERPOLATION_GEOMETRY) THEN
  583. ! MeshElementsConcThree=MeshElementsGeometry
  584. ! ! ! MESH_COMPONENT_NUMBER_CONC_THREE=MESH_COMPONENT_NUMBER_GEOMETRY
  585. ! ! ! ELSE IF(BASIS_XI_INTERPOLATION_CONC_THREE==BASIS_XI_INTERPOLATION_CONC_ONE) THEN
  586. ! ! ! MeshElementsConcThree=MeshElementsConcOne
  587. ! ! ! MESH_COMPONENT_NUMBER_CONC_THREE=MESH_COMPONENT_NUMBER_CONC_ONE
  588. ! ! ! ELSE
  589. ! ! ! MESH_COMPONENT_NUMBER_CONC_TWO=MESH_COMPONENT_NUMBER_CONC_ONE+1
  590. ! ! CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_CONC_THREE,BasisConcThree,MeshElementsConcThree,Err)
  591. ! ! ! DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
  592. ! ! ! CALL CMISSMeshElements_NodesSet(MeshElementsConcTwo,ELEMENT_NUMBER,CM%P(ELEMENT_NUMBER, &
  593. ! ! ! & 1:NUMBER_OF_ELEMENT_NODES_CONC_TWO),Err)
  594. ! ! ! ENDDO
  595. ! ! CALL CMISSMeshElements_CreateFinish(MeshElementsConcThree,Err)
  596. ! ! ! ENDIF
  597. !
  598. ! !Finish the creation of the mesh
  599. ! CALL CMISSMesh_CreateFinish(Mesh,Err)
  600. !Start the creation of a generated mesh in the region
  601. CALL CMISSGeneratedMesh_Initialise(GeneratedMesh,Err)
  602. CALL CMISSGeneratedMesh_CreateStart(GeneratedMeshUserNumber,Region,GeneratedMesh,Err)
  603. !Set up a regular x*y*z mesh
  604. CALL CMISSGeneratedMesh_TypeSet(GeneratedMesh,CMISS_GENERATED_MESH_REGULAR_MESH_TYPE,Err)
  605. !Set the default basis
  606. CALL CMISSGeneratedMesh_BasisSet(GeneratedMesh,BasisGeometry,Err)
  607. !Define the mesh on the region
  608. IF(NUMBER_GLOBAL_Z_ELEMENTS==0) THEN
  609. CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/WIDTH,HEIGHT/),Err)
  610. CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS/),Err)
  611. ELSE
  612. CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/WIDTH,HEIGHT,LENGTH/),Err)
  613. CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS, &
  614. & NUMBER_GLOBAL_Z_ELEMENTS/),Err)
  615. ENDIF
  616. !Finish the creation of a generated mesh in the region
  617. CALL CMISSMesh_Initialise(Mesh,Err)
  618. CALL CMISSGeneratedMesh_CreateFinish(GeneratedMesh,MeshUserNumber,Mesh,Err)
  619. !
  620. !================================================================================================================================
  621. !
  622. !GEOMETRIC FIELD
  623. !Create a decomposition:
  624. !All mesh components share the same decomposition
  625. CALL CMISSDecomposition_Initialise(Decomposition,Err)
  626. CALL CMISSDecomposition_CreateStart(DecompositionUserNumber,Mesh,Decomposition,Err)
  627. !Set the decomposition to be a general decomposition with the specified number of domains
  628. CALL CMISSDecomposition_TypeSet(Decomposition,CMISS_DECOMPOSITION_CALCULATED_TYPE,Err)
  629. CALL CMISSDecomposition_NumberOfDomainsSet(Decomposition,NumberOfUserDomains,Err)
  630. !Finish the decomposition
  631. CALL CMISSDecomposition_CreateFinish(Decomposition,Err)
  632. !Start to create a default (geometric) field on the region
  633. CALL CMISSField_Initialise(GeometricField,Err)
  634. CALL CMISSField_CreateStart(GeometricFieldUserNumber,Region,GeometricField,Err)
  635. !Set the field type
  636. CALL CMISSField_TypeSet(GeometricField,CMISS_FIELD_GEOMETRIC_TYPE,Err)
  637. !Set the decomposition to use
  638. CALL CMISSField_MeshDecompositionSet(GeometricField,Decomposition,Err)
  639. !Set the scaling to use
  640. CALL CMISSField_ScalingTypeSet(GeometricField,CMISS_FIELD_NO_SCALING,Err)
  641. DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  642. CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  643. & MESH_COMPONENT_NUMBER_GEOMETRY,Err)
  644. ENDDO
  645. !Finish creating the field
  646. CALL CMISSField_CreateFinish(GeometricField,Err)
  647. !Update the geometric field parameters
  648. CALL CMISSGeneratedMesh_GeometricParametersCalculate(GeneratedMesh,GeometricField,Err)
  649. !Update the geometric field parameters
  650. ! DO NODE_NUMBER=1,NUMBER_OF_NODES_GEOMETRY
  651. ! DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  652. ! VALUE=CM%N(NODE_NUMBER,COMPONENT_NUMBER)
  653. ! CALL CMISSField_ParameterSetUpdateNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  654. ! & CMISS_NO_GLOBAL_DERIV,NODE_NUMBER,COMPONENT_NUMBER,VALUE,Err)
  655. ! ENDDO
  656. ! ENDDO
  657. ! CALL CMISSField_ParameterSetUpdateStart(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,Err)
  658. ! CALL CMISSField_ParameterSetUpdateFinish(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,Err)
  659. !
  660. !================================================================================================================================
  661. !
  662. !EQUATIONS SETS - USING NEW ARGUMENTS TO ALLOW FOR MULTI-COMPARTMENT MODELS
  663. DO icompartment = 1,Ncompartments
  664. EquationsSetFieldUserNumberDiffusion = 100_CMISSIntg+icompartment
  665. EquationsSetUserNumberDiffusion = 200_CMISSIntg+icompartment
  666. CALL CMISSField_Initialise(EquationsSetFieldDiffusion(icompartment),Err)
  667. CALL CMISSEquationsSet_Initialise(EquationsSetDiffusion(icompartment),Err)
  668. CALL CMISSEquationsSet_CreateStart(EquationsSetUserNumberDiffusion,Region,GeometricField,&
  669. & CMISS_EQUATIONS_SET_CLASSICAL_FIELD_CLASS, &
  670. & CMISS_EQUATIONS_SET_DIFFUSION_EQUATION_TYPE,CMISS_EQUATIONS_SET_MULTI_COMP_TRANSPORT_DIFFUSION_SUBTYPE,&
  671. & EquationsSetFieldUserNumberDiffusion,EquationsSetFieldDiffusion(icompartment),EquationsSetDiffusion(icompartment),&
  672. & Err)
  673. !Finish creating the equations set
  674. CALL CMISSEquationsSet_CreateFinish(EquationsSetDiffusion(icompartment),Err)
  675. CALL CMISSField_ParameterSetUpdateConstant(RegionUserNumber,EquationsSetFielduserNumberDiffusion,CMISS_FIELD_U_VARIABLE_TYPE, &
  676. & CMISS_FIELD_VALUES_SET_TYPE,1,icompartment,Err)
  677. CALL CMISSField_ParameterSetUpdateConstant(EquationsSetFieldDiffusion(icompartment),CMISS_FIELD_U_VARIABLE_TYPE, &
  678. & CMISS_FIELD_VALUES_SET_TYPE,2,Ncompartments,Err)
  679. END DO
  680. !-------------------------------------------------------------------------------------
  681. ! DEPENDENT FIELD: Shared
  682. CALL CMISSField_Initialise(DependentField,Err)
  683. CALL CMISSField_CreateStart(DependentFieldUserNumber,Region,DependentField,Err)
  684. CALL CMISSField_TypeSet(DependentField,CMISS_FIELD_GENERAL_TYPE,Err)
  685. CALL CMISSField_MeshDecompositionSet(DependentField,Decomposition,Err)
  686. CALL CMISSField_GeometricFieldSet(DependentField,GeometricField,Err)
  687. CALL CMISSField_DependentTypeSet(DependentField,CMISS_FIELD_DEPENDENT_TYPE,Err)
  688. !Create 2N number of variables
  689. CALL CMISSField_NumberOfVariablesSet(DependentField,2*Ncompartments,Err)
  690. !create two variables for each compartment
  691. ALLOCATE(VariableTypes(2*Ncompartments))
  692. DO num_var=1,Ncompartments
  693. VariableTypes(2*num_var-1)=CMISS_FIELD_U_VARIABLE_TYPE+(CMISS_FIELD_NUMBER_OF_VARIABLE_SUBTYPES*(num_var-1))
  694. VariableTypes(2*num_var)=CMISS_FIELD_DELUDELN_VARIABLE_TYPE+(CMISS_FIELD_NUMBER_OF_VARIABLE_SUBTYPES*(num_var-1))
  695. ENDDO
  696. CALL CMISSField_VariableTypesSet(DependentField,VariableTypes,Err)
  697. !loop over the number of compartments
  698. DO icompartment=1,2*Ncompartments
  699. !set dimension type
  700. CALL CMISSField_DimensionSet(DependentField,VariableTypes(icompartment), &
  701. & CMISS_FIELD_SCALAR_DIMENSION_TYPE,Err)
  702. CALL CMISSField_NumberOfComponentsSet(DependentField,VariableTypes(icompartment),1,Err)
  703. CALL CMISSField_ComponentMeshComponentSet(DependentField,VariableTypes(icompartment),1, &
  704. & MESH_COMPONENT_NUMBER_CONC_ONE,Err)
  705. ENDDO
  706. CALL CMISSField_CreateFinish(DependentField,Err)
  707. DO icompartment = 1,Ncompartments
  708. CALL CMISSEquationsSet_DependentCreateStart(EquationsSetDiffusion(icompartment),DependentFieldUserNumber,&
  709. & DependentField,Err)
  710. CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetDiffusion(icompartment),Err)
  711. ENDDO
  712. ! !Set the mesh component to be used by the field components.
  713. ! ! DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  714. ! ! CALL CMISSField_ComponentMeshComponentSet(DependentFieldDiffusionTwo,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  715. ! ! & MESH_COMPONENT_NUMBER_CONC_TWO,Err)
  716. ! ! CALL CMISSField_ComponentMeshComponentSet(DependentFieldDiffusionTwo,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, &
  717. ! ! & MESH_COMPONENT_NUMBER_CONC_TWO,Err)
  718. ! ! ENDDO
  719. ! !Set the mesh component to be used by the field components.
  720. ! ! DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  721. ! ! CALL CMISSField_ComponentMeshComponentSet(DependentFieldDiffusionTwo,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  722. ! ! & MESH_COMPONENT_NUMBER_CONC_TWO,Err)
  723. ! ! CALL CMISSField_ComponentMeshComponentSet(DependentFieldDiffusionTwo,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, &
  724. ! ! & MESH_COMPONENT_NUMBER_CONC_TWO,Err)
  725. ! ! ENDDO
  726. ! !Set the mesh component to be used by the field components.
  727. ! ! DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
  728. ! ! CALL CMISSField_ComponentMeshComponentSet(DependentFieldDiffusionTwo,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, &
  729. ! ! & MESH_COMPONENT_NUMBER_CONC_TWO,Err)
  730. ! ! CALL CMISSField_ComponentMeshComponentSet(DependentFieldDiffusionTwo,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, &
  731. ! ! & MESH_COMPONENT_NUMBER_CONC_TWO,Err)
  732. ! ! ENDDO
  733. ! !Finish the equations set dependent field variables
  734. !-------------------------------------------------------------------------------------
  735. ! INITIALISE DEPENDENT FIELDS
  736. ! !Initialise dependent field (concentration one components)
  737. ! CALL CMISSField_ComponentValuesInitialise(DependentField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  738. ! & 1,1.0_CMISSDP,Err)
  739. ! !Initialise dependent field (concentration two components)
  740. ! CALL CMISSField_ComponentValuesInitialise(DependentField,CMISS_FIELD_U2_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  741. ! & 1,INITIAL_FIELD_DIFFUSION_TWO,Err)
  742. ! !Initialise dependent field (concentration three components)
  743. ! CALL CMISSField_ComponentValuesInitialise(DependentField,CMISS_FIELD_U3_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  744. ! & 1,INITIAL_FIELD_DIFFUSION_THREE,Err)
  745. !
  746. !================================================================================================================================
  747. !
  748. ALLOCATE(CouplingCoeffs(Ncompartments,Ncompartments))
  749. IF(Ncompartments==2)THEN
  750. ! CouplingCoeffs(1,1)=1.0E-01_CMISSDP
  751. ! ! CouplingCoeffs(1,2)=-1.0E-04_CMISSDP
  752. ! ! CouplingCoeffs(2,1)=-1.0E-04_CMISSDP
  753. ! CouplingCoeffs(1,2)=1.0E-01_CMISSDP
  754. ! CouplingCoeffs(2,1)=1.0E-01_CMISSDP
  755. ! CouplingCoeffs(2,2)=1.0E-01_CMISSDP
  756. CouplingCoeffs(1,1)=1.0E-01_CMISSDP
  757. CouplingCoeffs(1,2)=1.0E-01_CMISSDP
  758. CouplingCoeffs(2,1)=1.0E-01_CMISSDP
  759. CouplingCoeffs(2,2)=1.0E-01_CMISSDP
  760. ELSE IF(Ncompartments==3)THEN
  761. CouplingCoeffs(1,1)=1.0E-02_CMISSDP
  762. CouplingCoeffs(1,2)=1.0E-02_CMISSDP
  763. CouplingCoeffs(1,3)=0.0E-02_CMISSDP
  764. CouplingCoeffs(2,1)=1.0E-02_CMISSDP
  765. CouplingCoeffs(2,2)=2.0E-02_CMISSDP
  766. CouplingCoeffs(2,3)=1.0E-02_CMISSDP
  767. CouplingCoeffs(3,1)=0.0E-02_CMISSDP
  768. CouplingCoeffs(3,2)=1.0E-02_CMISSDP
  769. CouplingCoeffs(3,3)=1.0E-02_CMISSDP
  770. ELSE IF(Ncompartments==4)THEN
  771. CouplingCoeffs(1,1)=0.0E-02_CMISSDP
  772. CouplingCoeffs(1,2)=0.0E-02_CMISSDP
  773. CouplingCoeffs(1,3)=0.0E-02_CMISSDP
  774. CouplingCoeffs(1,4)=0.0E-02_CMISSDP
  775. CouplingCoeffs(2,1)=0.0E-02_CMISSDP
  776. CouplingCoeffs(2,2)=0.0E-02_CMISSDP
  777. CouplingCoeffs(2,3)=0.0E-02_CMISSDP
  778. CouplingCoeffs(2,4)=0.0E-02_CMISSDP
  779. CouplingCoeffs(3,1)=0.0E-02_CMISSDP
  780. CouplingCoeffs(3,2)=0.0E-02_CMISSDP
  781. CouplingCoeffs(3,3)=0.0E-02_CMISSDP
  782. CouplingCoeffs(3,4)=0.0E-02_CMISSDP
  783. CouplingCoeffs(4,1)=0.0E-02_CMISSDP
  784. CouplingCoeffs(4,2)=0.0E-02_CMISSDP
  785. CouplingCoeffs(4,3)=0.0E-02_CMISSDP
  786. CouplingCoeffs(4,4)=0.0E-02_CMISSDP
  787. ELSE IF(Ncompartments==5)THEN
  788. CouplingCoeffs(1,1)=0.0E-02_CMISSDP
  789. CouplingCoeffs(1,2)=0.0E-02_CMISSDP
  790. CouplingCoeffs(1,3)=0.0E-02_CMISSDP
  791. CouplingCoeffs(1,4)=0.0E-02_CMISSDP
  792. CouplingCoeffs(1,5)=0.0E-02_CMISSDP
  793. CouplingCoeffs(2,1)=0.0E-02_CMISSDP
  794. CouplingCoeffs(2,2)=0.0E-02_CMISSDP
  795. CouplingCoeffs(2,3)=0.0E-02_CMISSDP
  796. CouplingCoeffs(2,4)=0.0E-02_CMISSDP
  797. CouplingCoeffs(2,5)=0.0E-02_CMISSDP
  798. CouplingCoeffs(3,1)=0.0E-02_CMISSDP
  799. CouplingCoeffs(3,2)=0.0E-02_CMISSDP
  800. CouplingCoeffs(3,3)=0.0E-02_CMISSDP
  801. CouplingCoeffs(3,4)=0.0E-02_CMISSDP
  802. CouplingCoeffs(3,5)=0.0E-02_CMISSDP
  803. CouplingCoeffs(4,1)=0.0E-02_CMISSDP
  804. CouplingCoeffs(4,2)=0.0E-02_CMISSDP
  805. CouplingCoeffs(4,3)=0.0E-02_CMISSDP
  806. CouplingCoeffs(4,4)=0.0E-02_CMISSDP
  807. CouplingCoeffs(4,5)=0.0E-02_CMISSDP
  808. CouplingCoeffs(5,1)=0.0E-02_CMISSDP
  809. CouplingCoeffs(5,2)=0.0E-02_CMISSDP
  810. CouplingCoeffs(5,3)=0.0E-02_CMISSDP
  811. CouplingCoeffs(5,4)=0.0E-02_CMISSDP
  812. CouplingCoeffs(5,5)=0.0E-02_CMISSDP
  813. ELSE
  814. write(*,*) "Can't initialise coupling coefficients array."
  815. ENDIF
  816. !MATERIALS FIELDS - create the materials field
  817. !Auto-created field contains a U variable type to store the diffusion coefficient(s)
  818. !It also contains a V variable type to store the coupling coefficients
  819. DO icompartment = 1,Ncompartments
  820. MaterialsFieldUserNumberDiffusion = 400+icompartment
  821. CALL CMISSField_Initialise(MaterialsFieldDiffusion(icompartment),Err)
  822. CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetDiffusion(icompartment),MaterialsFieldUserNumberDiffusion,&
  823. & MaterialsFieldDiffusion(icompartment),Err)
  824. CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetDiffusion(icompartment),Err)
  825. END DO
  826. !Initialise the coupling coefficients
  827. !Need to devise a neater way of specifying these components - e.g. specify only the upper diagonal components, and then automatically fill out the rest
  828. DO icompartment = 1, Ncompartments
  829. DO COMPONENT_NUMBER=1, Ncompartments
  830. ! CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDiffusion(icompartment),CMISS_FIELD_V_VARIABLE_TYPE, &
  831. ! & CMISS_FIELD_VALUES_SET_TYPE,COMPONENT_NUMBER,CouplingCoeffs(icompartment,COMPONENT_NUMBER),Err)
  832. CALL CMISSField_ParameterSetUpdateConstant(MaterialsFieldDiffusion(icompartment),CMISS_FIELD_V_VARIABLE_TYPE, &
  833. & CMISS_FIELD_VALUES_SET_TYPE,COMPONENT_NUMBER,CouplingCoeffs(icompartment,COMPONENT_NUMBER),Err)
  834. END DO
  835. END DO
  836. ! CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDiffusionOne,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  837. ! & MaterialsFieldUserNumberDiffusionOne,POROSITY_PARAM_MAT_PROPERTIES,Err)
  838. ! CALL CMISSField_ComponentValuesInitialise(MaterialsFieldDiffusionTwo,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, &
  839. ! & MaterialsFieldUserNumberMatPropertiesPermOverVis,PERM_OVER_VIS_PARAM_MAT_PROPERTIES,Err)
  840. !
  841. !================================================================================================================================
  842. !
  843. !SOURCE FIELDS
  844. !create the equations set source field variables for both equations sets
  845. DO icompartment = 1,Ncompartments
  846. SourceFieldUserNumberDiffusion=700_CMISSIntg+icompartment
  847. CALL CMISSField_Initialise(SourceFieldDiffusion(icompartment),Err)
  848. CALL CMISSEquationsSet_SourceCreateStart(EquationsSetDiffusion(icompartment),SourceFieldUserNumberDiffusion,&
  849. & SourceFieldDiffusion(icompartment),Err)
  850. CALL CMISSEquationsSet_SourceCreateFinish(EquationsSetDiffusion(icompartment),Err)
  851. END DO
  852. !Create the equations set analytic field variables
  853. DO icompartment = 1,Ncompartments
  854. AnalyticFieldUserNumberDiffusion=900_CMISSIntg+icompartment
  855. CALL CMISSField_Initialise(AnalyticFieldDiffusion(icompartment),Err)
  856. CALL CMISSEquationsSet_AnalyticCreateStart(EquationsSetDiffusion(icompartment), &
  857. & CMISS_EQUATIONS_SET_MULTI_COMP_DIFFUSION_TWO_COMP_TWO_DIM,&
  858. & AnalyticFieldUserNumberDiffusion,AnalyticFieldDiffusion(icompartment),Err)
  859. !Finish the equations set analytic field variables
  860. CALL CMISSEquationsSet_AnalyticCreateFinish(EquationsSetDiffusion(icompartment),Err)
  861. END DO
  862. !
  863. !================================================================================================================================
  864. !
  865. !EQUATIONS
  866. DO icompartment=1,Ncompartments
  867. CALL CMISSEquations_Initialise(EquationsDiffusion(icompartment),Err)
  868. CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetDiffusion(icompartment),EquationsDiffusion(icompartment),Err)
  869. CALL CMISSEquations_SparsityTypeSet(EquationsDiffusion(icompartment),CMISS_EQUATIONS_SPARSE_MATRICES,Err)
  870. CALL CMISSEquations_OutputTypeSet(EquationsDiffusion(icompartment),EQUATIONS_DIFFUSION_OUTPUT,Err)
  871. CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetDiffusion(icompartment),Err)
  872. ENDDO
  873. !
  874. !================================================================================================================================
  875. !
  876. !PROBLEMS
  877. !Start the creation of a problem.
  878. CALL CMISSProblem_Initialise(Problem,Err)
  879. CALL CMISSControlLoop_Initialise(ControlLoop,Err)
  880. CALL CMISSProblem_CreateStart(ProblemUserNumber,Problem,Err)
  881. !Set the problem to be a coupled diffusion-diffusion problem
  882. CALL CMISSProblem_SpecificationSet(Problem,CMISS_PROBLEM_MULTI_PHYSICS_CLASS,CMISS_PROBLEM_MULTI_COMPARTMENT_TRANSPORT_TYPE, &
  883. & CMISS_PROBLEM_STANDARD_MULTI_COMPARTMENT_TRANSPORT_SUBTYPE,Err)
  884. !Finish the creation of a problem.
  885. CALL CMISSProblem_CreateFinish(Problem,Err)
  886. !Start the creation of the problem control loop
  887. CALL CMISSProblem_ControlLoopCreateStart(Problem,Err)
  888. !Get the control loop
  889. CALL CMISSProblem_ControlLoopGet(Problem,CMISS_CONTROL_LOOP_NODE,ControlLoop,Err)
  890. !Set the times
  891. CALL CMISSControlLoop_TimesSet(ControlLoop,LINEAR_SOLVER_DIFFUSION_START_TIME,&
  892. & LINEAR_SOLVER_DIFFUSION_STOP_TIME,LINEAR_SOLVER_DIFFUSION_TIME_INCREMENT,Err)
  893. !Set the output timing
  894. CALL CMISSControlLoop_TimeOutputSet(ControlLoop,LINEAR_SOLVER_DIFFUSION_OUTPUT_FREQUENCY,Err)
  895. !Finish creating the problem control loop
  896. CALL CMISSProblem_ControlLoopCreateFinish(Problem,Err)
  897. !
  898. !================================================================================================================================
  899. !
  900. !SOLVERS
  901. !Start the creation of the problem solvers
  902. CALL CMISSSolver_Initialise(SolverDiffusion,Err)
  903. CALL CMISSSolver_Initialise(LinearSolverDiffusion,Err)
  904. CALL CMISSProblem_SolversCreateStart(Problem,Err)
  905. CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,SolverDiffusionUserNumber,SolverDiffusion,Err)
  906. !Set the output type
  907. CALL CMISSSolver_OutputTypeSet(SolverDiffusion,LINEAR_SOLVER_DIFFUSION_OUTPUT_TYPE,Err)
  908. !Set the solver settings
  909. ! IF(LINEAR_SOLVER_DIFFUSION_ONE_DIRECT_FLAG) THEN
  910. ! CALL CMISSSolver_LinearTypeSet(LinearSolverDiffusionOne,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err)
  911. ! CALL CMISSSolver_LibraryTypeSet(LinearSolverDiffusionOne,CMISS_SOLVER_MUMPS_LIBRARY,Err)
  912. ! ELSE
  913. ! CALL CMISSSolver_LinearTypeSet(LinearSolverDiffusionOne,CMISS_SOLVER_LINEAR_ITERATIVE_SOLVE_TYPE,Err)
  914. CALL CMISSSolver_DynamicLinearSolverGet(SolverDiffusion,LinearSolverDiffusion,Err)
  915. CALL CMISSSolver_LinearIterativeMaximumIterationsSet(LinearSolverDiffusion,MAXIMUM_ITERATIONS,Err)
  916. ! CALL CMISSSolver_LinearIterativeDivergenceToleranceSet(LinearSolverDiffusionOne,DIVERGENCE_TOLERANCE,Err)
  917. ! CALL CMISSSolver_LinearIterativeRelativeToleranceSet(LinearSolverDiffusionOne,RELATIVE_TOLERANCE,Err)
  918. ! CALL CMISSSolver_LinearIterativeAbsoluteToleranceSet(LinearSolverDiffusionOne,ABSOLUTE_TOLERANCE,Err)
  919. ! CALL CMISSSolver_LinearIterativeGMRESRestartSet(LinearSolverDiffusionOne,RESTART_VALUE,Err)
  920. ! ENDIF
  921. !Finish the creation of the problem solver
  922. CALL CMISSProblem_SolversCreateFinish(Problem,Err)
  923. !
  924. !================================================================================================================================
  925. !
  926. !SOLVER EQUATIONS
  927. !Start the creation of the problem solver equations
  928. CALL CMISSSolver_Initialise(SolverDiffusion,Err)
  929. CALL CMISSSolverEquations_Initialise(SolverEquationsDiffusion,Err)
  930. CALL CMISSProblem_SolverEquationsCreateStart(Problem,Err)
  931. !
  932. CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,SolverDiffusionUserNumber,SolverDiffusion,Err)
  933. CALL CMISSSolver_SolverEquationsGet(SolverDiffusion,SolverEquationsDiffusion,Err)
  934. CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsDiffusion,CMISS_SOLVER_SPARSE_MATRICES,Err)
  935. DO icompartment=1,Ncompartments
  936. CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsDiffusion,EquationsSetDiffusion(icompartment),&
  937. & EquationsSetIndex,Err)
  938. ENDDO
  939. !Finish the creation of the problem solver equations
  940. CALL CMISSProblem_SolverEquationsCreateFinish(Problem,Err)
  941. !
  942. !================================================================================================================================
  943. !
  944. !BOUNDARY CONDITIONS
  945. CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsDiffusion,BoundaryConditionsDiffusion,Err)
  946. CALL CMISSSolverEquations_BoundaryConditionsAnalytic(SolverEquationsDiffusion,Err)
  947. CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquationsDiffusion,Err)
  948. !
  949. !================================================================================================================================
  950. !
  951. !RUN SOLVERS
  952. !Turn off PETSc error handling
  953. !CALL PETSC_ERRORHANDLING_SET_ON(ERR,ERROR,*999)
  954. !Solve the problem
  955. CALL CMISSProblem_Solve(Problem,Err)
  956. !
  957. !================================================================================================================================
  958. !
  959. !OUTPUT
  960. Call CMISSAnalyticAnalysisOutput(DependentField,"MultiCompDiffusionAnalytics_2Comp_2D_x40_y40_L_T1",Err)
  961. EXPORT_FIELD_IO=.TRUE.
  962. IF(EXPORT_FIELD_IO) THEN
  963. WRITE(*,'(A)') "Exporting fields..."
  964. ! CALL CMISSFields_Initialise(Fields,Err)
  965. ! CALL CMISSFields_Create(Region,Fields,Err)
  966. ! CALL CMISSFields_NodesExport(Fields,"MonolithicMultiCompDiffusionTest","FORTRAN",Err)
  967. ! CALL CMISSFields_ElementsExport(Fields,"MonolithicMultiCompDiffusionTest","FORTRAN",Err)
  968. ! CALL CMISSFields_Finalise(Fields,Err)
  969. WRITE(*,'(A)') "Field exported!"
  970. CALL FieldmlOutput_InitializeInfo( Region, Mesh, 3, outputDirectory, basename, fieldmlInfo, err )
  971. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".geometric", region, mesh, GeometricField, &
  972. & CMISS_FIELD_U_VARIABLE_TYPE, err )
  973. ! CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".dependent", region, mesh, DependentField, err )
  974. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".dependent_U", region, mesh, DependentField, &
  975. & CMISS_FIELD_U_VARIABLE_TYPE, err )
  976. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".dependent_V", region, mesh, DependentField, &
  977. & CMISS_FIELD_V_VARIABLE_TYPE, err )
  978. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".equations_set_field_1", region, mesh, EquationsSetFieldDiffusion(1), &
  979. & CMISS_FIELD_U_VARIABLE_TYPE, err )
  980. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".equations_set_field_2", region, mesh, EquationsSetFieldDiffusion(2), &
  981. & CMISS_FIELD_U_VARIABLE_TYPE, err )
  982. ! CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".equations_set_field_3", region, mesh, EquationsSetFieldDiffusion(3), &
  983. ! & CMISS_FIELD_U_VARIABLE_TYPE, err )
  984. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".source_1", region,mesh,SourceFieldDiffusion(1), &
  985. & CMISS_FIELD_U_VARIABLE_TYPE,&
  986. & err )
  987. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".source_2", region,mesh,SourceFieldDiffusion(2), &
  988. & CMISS_FIELD_U_VARIABLE_TYPE,&
  989. & err )
  990. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".materials_1", region, mesh, MaterialsFieldDiffusion(1), &
  991. & CMISS_FIELD_U_VARIABLE_TYPE,err )
  992. CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".materials_2", region, mesh, MaterialsFieldDiffusion(2), &
  993. & CMISS_FIELD_U_VARIABLE_TYPE,err )
  994. !CALL FieldmlOutput_AddField( fieldmlInfo, baseName//".analytic", region, mesh, AnalyticField, err )
  995. CALL FieldmlOutput_Write( fieldmlInfo, outputFilename, err )
  996. CALL FieldmlUtil_FinalizeInfo( fieldmlInfo )
  997. ENDIF
  998. !Finialise CMISS
  999. ! CALL CMISSFinalise(Err)
  1000. WRITE(*,'(A)') "Program successfully completed."
  1001. STOP
  1002. END PROGRAM MONOLITHICSCHEMETESTFIELDMLEXAMPLE