PageRenderTime 105ms CodeModel.GetById 13ms app.highlight 84ms RepoModel.GetById 1ms app.codeStats 1ms

/FluidMechanics/NavierStokes/temp_tests/DynamicTest_theta_0.5/src/DynamicTest_theta_0.5Example.f90

http://github.com/xyan075/examples
FORTRAN Modern | 839 lines | 467 code | 110 blank | 262 comment | 0 complexity | 3c7a4dd0999b45c5fd0b072db7b19a2f MD5 | raw file
  1!> \file
  2!> \author Sebastian Krittian
  3!> \brief This is an example program to solve an analytic static Navier-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
 42!> \example FluidMechanics/NavierStokes/Static/src/StaticExample.f90
 43!! Example program to solve an analytic static Navier-Stokes equation using OpenCMISS calls.
 44!! \par Latest Builds:
 45!! \li <a href='http://autotest.bioeng.auckland.ac.nz/opencmiss-build/logs_x86_64-linux/FluidMechanics/NavierStokes/Static/build-intel'>Linux Intel Build</a>
 46!! \li <a href='http://autotest.bioeng.auckland.ac.nz/opencmiss-build/logs_x86_64-linux/FluidMechanics/NavierStokes/Static/build-intel'>Linux GNU Build</a>
 47!!
 48!<
 49
 50!> Main program
 51
 52PROGRAM ANALYTICNAVIERSTOKESEXAMPLE
 53
 54  !
 55  !================================================================================================================================
 56  !
 57
 58  !PROGRAM LIBRARIES
 59
 60  USE OPENCMISS
 61  USE FLUID_MECHANICS_IO_ROUTINES
 62  USE MPI
 63
 64#ifdef WIN32
 65  USE IFQWINCMISS
 66#endif
 67
 68  !
 69  !================================================================================================================================
 70  !
 71
 72  !PROGRAM VARIABLES AND TYPES
 73
 74  IMPLICIT NONE
 75
 76  INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldUserNumber=1337
 77  TYPE(CMISSFieldType) :: EquationsSetField
 78
 79
 80  !Test program parameters
 81
 82  INTEGER(CMISSIntg), PARAMETER :: CoordinateSystemUserNumber=1
 83  INTEGER(CMISSIntg), PARAMETER :: RegionUserNumber=2
 84  INTEGER(CMISSIntg), PARAMETER :: MeshUserNumber=3
 85  INTEGER(CMISSIntg), PARAMETER :: DecompositionUserNumber=4
 86  INTEGER(CMISSIntg), PARAMETER :: GeometricFieldUserNumber=5
 87  INTEGER(CMISSIntg), PARAMETER :: DependentFieldUserNumberNavierStokes=6
 88  INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberNavierStokes=7
 89  INTEGER(CMISSIntg), PARAMETER :: IndependentFieldUserNumberNavierStokes=8
 90  INTEGER(CMISSIntg), PARAMETER :: AnalyticFieldUserNumberNavierStokes=9
 91  INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumberNavierStokes=10
 92  INTEGER(CMISSIntg), PARAMETER :: ProblemUserNumber=11
 93
 94  INTEGER(CMISSIntg), PARAMETER :: DomainUserNumber=1
 95  INTEGER(CMISSIntg), PARAMETER :: SolverNavierStokesUserNumber=1
 96  INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberNavierStokesMu=1
 97  INTEGER(CMISSIntg), PARAMETER :: MaterialsFieldUserNumberNavierStokesRho=2
 98
 99  !Program types
100
101  TYPE(EXPORT_CONTAINER):: CM
102
103  !Program variables
104
105  INTEGER(CMISSIntg) :: NUMBER_OF_DIMENSIONS
106  
107  INTEGER(CMISSIntg) :: BASIS_TYPE
108  INTEGER(CMISSIntg) :: BASIS_NUMBER_SPACE
109  INTEGER(CMISSIntg) :: BASIS_NUMBER_VELOCITY
110  INTEGER(CMISSIntg) :: BASIS_NUMBER_PRESSURE
111  INTEGER(CMISSIntg) :: BASIS_GAUSS_SPACE
112  INTEGER(CMISSIntg) :: BASIS_GAUSS_VELOCITY
113  INTEGER(CMISSIntg) :: BASIS_GAUSS_PRESSURE
114  INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_SPACE
115  INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_VELOCITY
116  INTEGER(CMISSIntg) :: BASIS_XI_INTERPOLATION_PRESSURE
117  INTEGER(CMISSIntg) :: MESH_NUMBER_OF_COMPONENTS
118  INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_SPACE
119  INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_VELOCITY
120  INTEGER(CMISSIntg) :: MESH_COMPONENT_NUMBER_PRESSURE
121  INTEGER(CMISSIntg) :: NUMBER_OF_NODES_SPACE
122  INTEGER(CMISSIntg) :: NUMBER_OF_NODES_VELOCITY
123  INTEGER(CMISSIntg) :: NUMBER_OF_NODES_PRESSURE
124  INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_SPACE
125  INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_VELOCITY
126  INTEGER(CMISSIntg) :: NUMBER_OF_ELEMENT_NODES_PRESSURE
127  INTEGER(CMISSIntg) :: TOTAL_NUMBER_OF_NODES
128  INTEGER(CMISSIntg) :: TOTAL_NUMBER_OF_ELEMENTS
129  INTEGER(CMISSIntg) :: MAXIMUM_ITERATIONS
130  INTEGER(CMISSIntg) :: RESTART_VALUE
131!   INTEGER(CMISSIntg) :: MPI_IERROR
132
133  INTEGER(CMISSIntg) :: EQUATIONS_NAVIER_STOKES_OUTPUT
134  INTEGER(CMISSIntg) :: COMPONENT_NUMBER
135  INTEGER(CMISSIntg) :: NODE_NUMBER
136  INTEGER(CMISSIntg) :: ELEMENT_NUMBER
137!   INTEGER(CMISSIntg) :: NODE_COUNTER
138!   INTEGER(CMISSIntg) :: CONDITION
139
140  INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_NAVIER_STOKES_OUTPUT_FREQUENCY
141  INTEGER(CMISSIntg) :: DYNAMIC_SOLVER_NAVIER_STOKES_OUTPUT_TYPE
142  INTEGER(CMISSIntg) :: NONLINEAR_SOLVER_NAVIER_STOKES_OUTPUT_TYPE
143  INTEGER(CMISSIntg) :: LINEAR_SOLVER_NAVIER_STOKES_OUTPUT_TYPE
144
145  INTEGER(CMISSIntg) :: ANALYTICAL_TYPE
146  INTEGER(CMISSIntg) :: INPUT_TYPE
147
148  REAL(CMISSDP) :: INITIAL_FIELD_NAVIER_STOKES(3)
149  REAL(CMISSDP) :: DIVERGENCE_TOLERANCE
150  REAL(CMISSDP) :: RELATIVE_TOLERANCE
151  REAL(CMISSDP) :: ABSOLUTE_TOLERANCE
152  REAL(CMISSDP) :: SOLVER_TOLERANCE
153  REAL(CMISSDP) :: LINESEARCH_ALPHA
154  REAL(CMISSDP) :: MAX_STEP
155  REAL(CMISSDP) :: VALUE
156  REAL(CMISSDP) :: MU_PARAM_NAVIER_STOKES
157  REAL(CMISSDP) :: RHO_PARAM_NAVIER_STOKES
158
159  REAL(CMISSDP) :: DYNAMIC_SOLVER_NAVIER_STOKES_START_TIME
160  REAL(CMISSDP) :: DYNAMIC_SOLVER_NAVIER_STOKES_STOP_TIME
161  REAL(CMISSDP) :: DYNAMIC_SOLVER_NAVIER_STOKES_THETA
162  REAL(CMISSDP) :: DYNAMIC_SOLVER_NAVIER_STOKES_TIME_INCREMENT
163
164  LOGICAL :: EXPORT_FIELD_IO
165  LOGICAL :: LINEAR_SOLVER_NAVIER_STOKES_DIRECT_FLAG
166
167  CHARACTER *15 BUFFER
168  CHARACTER *15 OUTPUT_STRING
169
170  !CMISS variables
171
172  !Regions
173  TYPE(CMISSRegionType) :: Region
174  TYPE(CMISSRegionType) :: WorldRegion
175  !Coordinate systems
176  TYPE(CMISSCoordinateSystemType) :: CoordinateSystem
177  TYPE(CMISSCoordinateSystemType) :: WorldCoordinateSystem
178  !Basis
179  TYPE(CMISSBasisType) :: BasisSpace
180  TYPE(CMISSBasisType) :: BasisVelocity
181  TYPE(CMISSBasisType) :: BasisPressure
182  !Nodes
183  TYPE(CMISSNodesType) :: Nodes
184  !Elements
185  TYPE(CMISSMeshElementsType) :: MeshElementsSpace
186  TYPE(CMISSMeshElementsType) :: MeshElementsVelocity
187  TYPE(CMISSMeshElementsType) :: MeshElementsPressure
188  !Meshes
189  TYPE(CMISSMeshType) :: Mesh
190  !Decompositions
191  TYPE(CMISSDecompositionType) :: Decomposition
192  !Fields
193  TYPE(CMISSFieldsType) :: Fields
194  !Field types
195  TYPE(CMISSFieldType) :: GeometricField
196  TYPE(CMISSFieldType) :: DependentFieldNavierStokes
197  TYPE(CMISSFieldType) :: MaterialsFieldNavierStokes
198  TYPE(CMISSFieldType) :: AnalyticFieldNavierStokes
199  !Boundary conditions
200   TYPE(CMISSBoundaryConditionsType) :: BoundaryConditionsNavierStokes
201  !Equations sets
202  TYPE(CMISSEquationsSetType) :: EquationsSetNavierStokes
203  !Equations
204  TYPE(CMISSEquationsType) :: EquationsNavierStokes
205  !Problems
206  TYPE(CMISSProblemType) :: Problem
207  !Control loops
208  TYPE(CMISSControlLoopType) :: ControlLoop
209  !Solvers
210  TYPE(CMISSSolverType) :: DynamicSolverNavierStokes
211  TYPE(CMISSSolverType) :: NonlinearSolverNavierStokes
212  TYPE(CMISSSolverType) :: LinearSolverNavierStokes
213  !Solver equations
214  TYPE(CMISSSolverEquationsType) :: SolverEquationsNavierStokes
215
216#ifdef WIN32
217  !Quickwin type
218  LOGICAL :: QUICKWIN_STATUS=.FALSE.
219  TYPE(WINDOWCONFIG) :: QUICKWIN_WINDOW_CONFIG
220#endif
221  
222  !Generic CMISS variables
223  
224  INTEGER(CMISSIntg) :: EquationsSetIndex
225  INTEGER(CMISSIntg) :: Err
226  
227#ifdef WIN32
228  !Initialise QuickWin
229  QUICKWIN_WINDOW_CONFIG%TITLE="General Output" !Window title
230  QUICKWIN_WINDOW_CONFIG%NUMTEXTROWS=-1 !Max possible number of rows
231  QUICKWIN_WINDOW_CONFIG%MODE=QWIN$SCROLLDOWN
232  !Set the window parameters
233  QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
234  !If attempt fails set with system estimated values
235  IF(.NOT.QUICKWIN_STATUS) QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG)
236#endif
237
238  !
239  !================================================================================================================================
240  !
241
242  !PROBLEM CONTROL PANEL
243
244  !Import cmHeart mesh information
245  CALL FLUID_MECHANICS_IO_READ_CMHEART(CM,Err)  
246  BASIS_NUMBER_SPACE=CM%ID_M
247  BASIS_NUMBER_VELOCITY=CM%ID_V
248  BASIS_NUMBER_PRESSURE=CM%ID_P
249  NUMBER_OF_DIMENSIONS=CM%D
250  BASIS_TYPE=CM%IT_T
251  BASIS_XI_INTERPOLATION_SPACE=CM%IT_M
252  BASIS_XI_INTERPOLATION_VELOCITY=CM%IT_V
253  BASIS_XI_INTERPOLATION_PRESSURE=CM%IT_P
254  NUMBER_OF_NODES_SPACE=CM%N_M
255  NUMBER_OF_NODES_VELOCITY=CM%N_V
256  NUMBER_OF_NODES_PRESSURE=CM%N_P
257  TOTAL_NUMBER_OF_NODES=CM%N_T
258  TOTAL_NUMBER_OF_ELEMENTS=CM%E_T
259  NUMBER_OF_ELEMENT_NODES_SPACE=CM%EN_M
260  NUMBER_OF_ELEMENT_NODES_VELOCITY=CM%EN_V
261  NUMBER_OF_ELEMENT_NODES_PRESSURE=CM%EN_P
262  !Set initial values
263  INITIAL_FIELD_NAVIER_STOKES(1)=0.0_CMISSDP
264  INITIAL_FIELD_NAVIER_STOKES(2)=0.0_CMISSDP
265  INITIAL_FIELD_NAVIER_STOKES(3)=0.0_CMISSDP
266  !Set material parameters
267  MU_PARAM_NAVIER_STOKES=1.0_CMISSDP
268  RHO_PARAM_NAVIER_STOKES=1.0_CMISSDP 
269  !Set interpolation parameters
270  BASIS_GAUSS_SPACE=4
271  BASIS_GAUSS_VELOCITY=4
272  BASIS_GAUSS_PRESSURE=4
273  !Set output parameter
274  !(NoOutput/ProgressOutput/TimingOutput/SolverOutput/SolverMatrixOutput)
275  DYNAMIC_SOLVER_NAVIER_STOKES_OUTPUT_TYPE=CMISS_SOLVER_NO_OUTPUT
276  NONLINEAR_SOLVER_NAVIER_STOKES_OUTPUT_TYPE=CMISS_SOLVER_PROGRESS_OUTPUT
277  LINEAR_SOLVER_NAVIER_STOKES_OUTPUT_TYPE=CMISS_SOLVER_NO_OUTPUT
278  !(NoOutput/TimingOutput/MatrixOutput/ElementOutput)
279  EQUATIONS_NAVIER_STOKES_OUTPUT=CMISS_EQUATIONS_NO_OUTPUT
280  !Set time parameter
281  DYNAMIC_SOLVER_NAVIER_STOKES_START_TIME=0.0_CMISSDP
282  DYNAMIC_SOLVER_NAVIER_STOKES_STOP_TIME=1.001_CMISSDP 
283  DYNAMIC_SOLVER_NAVIER_STOKES_TIME_INCREMENT=0.100_CMISSDP
284  DYNAMIC_SOLVER_NAVIER_STOKES_THETA=1.0_CMISSDP/2.0_CMISSDP
285  !Set result output parameter
286  DYNAMIC_SOLVER_NAVIER_STOKES_OUTPUT_FREQUENCY=1
287  !Set solver parameters
288  LINEAR_SOLVER_NAVIER_STOKES_DIRECT_FLAG=.TRUE.
289  RELATIVE_TOLERANCE=1.0E-8_CMISSDP !default: 1.0E-05_CMISSDP
290  ABSOLUTE_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-10_CMISSDP
291  SOLVER_TOLERANCE=1.0E-10_CMISSDP !default: 1.0E-10_CMISSDP
292  DIVERGENCE_TOLERANCE=1.0E20 !default: 1.0E5
293  MAXIMUM_ITERATIONS=100000 !default: 100000
294  RESTART_VALUE=3000 !default: 30
295  LINESEARCH_ALPHA=1.0
296  MAX_STEP=0.00000001
297
298  !
299  !================================================================================================================================
300  !
301
302  !Different analytical cases
303  ! 1=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_1
304  ! 2=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_2
305  ! 3=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_3
306  ! 4=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_THREE_DIM_1
307  ! 5=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_THREE_DIM_2
308  ! 6=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_THREE_DIM_3
309
310
311! ! !   WRITE(*,*)'1=POLYNOM, 2=EXP, 3=COS/SIN:'
312! ! ! !   READ(*,*) 
313! ! ! 
314! ! !   IF(COMMAND_ARGUMENT_COUNT()==2) THEN
315! ! !     CALL GET_COMMAND_ARGUMENT(1,BUFFER)
316! ! !     READ(BUFFER,*) INPUT_TYPE
317! ! !     CALL GET_COMMAND_ARGUMENT(2,BUFFER)
318! ! !     READ(BUFFER,*) OUTPUT_STRING
319! ! !   ELSE
320! ! !     !TODO more detailed error message
321! ! !     WRITE(*,*)'INPUT ERROR!!!'
322! ! !   ENDIF
323
324  INPUT_TYPE=5
325  OUTPUT_STRING="TaylorGreen"
326
327
328  IF(INPUT_TYPE==1.AND.NUMBER_OF_DIMENSIONS==2) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_1
329  IF(INPUT_TYPE==2.AND.NUMBER_OF_DIMENSIONS==2) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_2
330  IF(INPUT_TYPE==3.AND.NUMBER_OF_DIMENSIONS==2) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_3
331  IF(INPUT_TYPE==4.AND.NUMBER_OF_DIMENSIONS==2) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_4
332  IF(INPUT_TYPE==5.AND.NUMBER_OF_DIMENSIONS==2) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TWO_DIM_5
333  IF(INPUT_TYPE==1.AND.NUMBER_OF_DIMENSIONS==3) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_THREE_DIM_1
334  IF(INPUT_TYPE==2.AND.NUMBER_OF_DIMENSIONS==3) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_THREE_DIM_2
335  IF(INPUT_TYPE==3.AND.NUMBER_OF_DIMENSIONS==3) ANALYTICAL_TYPE=CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_THREE_DIM_3
336
337  !
338  !================================================================================================================================
339  !
340
341  !INITIALISE OPENCMISS
342
343  CALL CMISSInitialise(WorldCoordinateSystem,WorldRegion,Err)
344
345  !
346  !================================================================================================================================
347  !
348
349  !COORDINATE SYSTEM
350
351  !Start the creation of a new RC coordinate system
352  CALL CMISSCoordinateSystem_Initialise(CoordinateSystem,Err)
353  CALL CMISSCoordinateSystem_CreateStart(CoordinateSystemUserNumber,CoordinateSystem,Err)
354  !Set the coordinate system dimension
355  CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,NUMBER_OF_DIMENSIONS,Err)
356  !Finish the creation of the coordinate system
357  CALL CMISSCoordinateSystem_CreateFinish(CoordinateSystem,Err)
358
359  !
360  !================================================================================================================================
361  !
362
363  !REGION
364
365  !Start the creation of a new region
366  CALL CMISSRegion_Initialise(Region,Err)
367  CALL CMISSRegion_CreateStart(RegionUserNumber,WorldRegion,Region,Err)
368  !Set the regions coordinate system as defined above
369  CALL CMISSRegion_CoordinateSystemSet(Region,CoordinateSystem,Err)
370  !Finish the creation of the region
371  CALL CMISSRegion_CreateFinish(Region,Err)
372
373  !
374  !================================================================================================================================
375  !
376
377  !BASES
378
379  !Start the creation of new bases
380  MESH_NUMBER_OF_COMPONENTS=1
381  CALL CMISSBasis_Initialise(BasisSpace,Err)
382  CALL CMISSBasis_CreateStart(BASIS_NUMBER_SPACE,BasisSpace,Err)
383  !Set the basis type (Lagrange/Simplex)
384  CALL CMISSBasis_TypeSet(BasisSpace,BASIS_TYPE,Err)
385  !Set the basis xi number
386  CALL CMISSBasis_NumberOfXiSet(BasisSpace,NUMBER_OF_DIMENSIONS,Err)
387  !Set the basis xi interpolation and number of Gauss points
388  IF(NUMBER_OF_DIMENSIONS==2) THEN
389    CALL CMISSBasis_InterpolationXiSet(BasisSpace,(/BASIS_XI_INTERPOLATION_SPACE,BASIS_XI_INTERPOLATION_SPACE/),Err)
390    IF(BASIS_TYPE/=CMISS_BASIS_SIMPLEX_TYPE) THEN
391      CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisSpace,(/BASIS_GAUSS_SPACE,BASIS_GAUSS_SPACE/),Err)
392    ELSE
393      CALL CMISSBasis_QuadratureOrderSet(BasisSpace,BASIS_GAUSS_SPACE+1,Err)
394    ENDIF
395  ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
396    CALL CMISSBasis_InterpolationXiSet(BasisSpace,(/BASIS_XI_INTERPOLATION_SPACE,BASIS_XI_INTERPOLATION_SPACE, & 
397      & BASIS_XI_INTERPOLATION_SPACE/),Err)                         
398    IF(BASIS_TYPE/=CMISS_BASIS_SIMPLEX_TYPE) THEN
399      CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisSpace,(/BASIS_GAUSS_SPACE,BASIS_GAUSS_SPACE,BASIS_GAUSS_SPACE/), & 
400        & Err)
401    ELSE
402      CALL CMISSBasis_QuadratureOrderSet(BasisSpace,BASIS_GAUSS_SPACE+1,Err)
403    ENDIF
404  ENDIF
405  !Finish the creation of the basis
406  CALL CMISSBasis_CreateFinish(BasisSpace,Err)
407  !Start the creation of another basis
408  IF(BASIS_XI_INTERPOLATION_VELOCITY==BASIS_XI_INTERPOLATION_SPACE) THEN
409    BasisVelocity=BasisSpace
410  ELSE
411    MESH_NUMBER_OF_COMPONENTS=MESH_NUMBER_OF_COMPONENTS+1
412    !Initialise a new velocity basis
413    CALL CMISSBasis_Initialise(BasisVelocity,Err)
414    !Start the creation of a basis
415    CALL CMISSBasis_CreateStart(BASIS_NUMBER_VELOCITY,BasisVelocity,Err)
416    !Set the basis type (Lagrange/Simplex)LinearSolverNavierStokes
417    CALL CMISSBasis_TypeSet(BasisVelocity,BASIS_TYPE,Err)
418    !Set the basis xi number
419    CALL CMISSBasis_NumberOfXiSet(BasisVelocity,NUMBER_OF_DIMENSIONS,Err)
420    !Set the basis xi interpolation and number of Gauss points
421    IF(NUMBER_OF_DIMENSIONS==2) THEN
422      CALL CMISSBasis_InterpolationXiSet(BasisVelocity,(/BASIS_XI_INTERPOLATION_VELOCITY,BASIS_XI_INTERPOLATION_VELOCITY/),Err)
423      IF(BASIS_TYPE/=CMISS_BASIS_SIMPLEX_TYPE) THEN 
424        CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisVelocity,(/BASIS_GAUSS_VELOCITY,BASIS_GAUSS_VELOCITY/),Err)
425      ELSE
426        CALL CMISSBasis_QuadratureOrderSet(BasisVelocity,BASIS_GAUSS_VELOCITY+1,Err)
427      ENDIF
428    ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
429      CALL CMISSBasis_InterpolationXiSet(BasisVelocity,(/BASIS_XI_INTERPOLATION_VELOCITY,BASIS_XI_INTERPOLATION_VELOCITY, & 
430        & BASIS_XI_INTERPOLATION_VELOCITY/),Err)                         
431      IF(BASIS_TYPE/=CMISS_BASIS_SIMPLEX_TYPE) THEN
432        CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisVelocity,(/BASIS_GAUSS_VELOCITY,BASIS_GAUSS_VELOCITY, & 
433          & BASIS_GAUSS_VELOCITY/),Err)
434      ELSE
435        CALL CMISSBasis_QuadratureOrderSet(BasisVelocity,BASIS_GAUSS_VELOCITY+1,Err)
436      ENDIF
437    ENDIF
438    !Finish the creation of the basis
439    CALL CMISSBasis_CreateFinish(BasisVelocity,Err)
440  ENDIF
441  !Start the creation of another basis
442  IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_SPACE) THEN
443    BasisPressure=BasisSpace
444  ELSE IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_VELOCITY) THEN
445    BasisPressure=BasisVelocity
446  ELSE
447    MESH_NUMBER_OF_COMPONENTS=MESH_NUMBER_OF_COMPONENTS+1
448    !Initialise a new pressure basis
449    CALL CMISSBasis_Initialise(BasisPressure,Err)
450    !Start the creation of a basis
451    CALL CMISSBasis_CreateStart(BASIS_NUMBER_PRESSURE,BasisPressure,Err)
452    !Set the basis type (Lagrange/Simplex)
453    CALL CMISSBasis_TypeSet(BasisPressure,BASIS_TYPE,Err)
454    !Set the basis xi number
455    CALL CMISSBasis_NumberOfXiSet(BasisPressure,NUMBER_OF_DIMENSIONS,Err)
456    !Set the basis xi interpolation and number of Gauss points
457    IF(NUMBER_OF_DIMENSIONS==2) THEN
458      CALL CMISSBasis_InterpolationXiSet(BasisPressure,(/BASIS_XI_INTERPOLATION_PRESSURE,BASIS_XI_INTERPOLATION_PRESSURE/),Err)
459      IF(BASIS_TYPE/=CMISS_BASIS_SIMPLEX_TYPE) THEN
460        CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisPressure,(/BASIS_GAUSS_PRESSURE,BASIS_GAUSS_PRESSURE/),Err)
461      ELSE
462        CALL CMISSBasis_QuadratureOrderSet(BasisPressure,BASIS_GAUSS_PRESSURE+1,Err)
463      ENDIF
464    ELSE IF(NUMBER_OF_DIMENSIONS==3) THEN
465      CALL CMISSBasis_InterpolationXiSet(BasisPressure,(/BASIS_XI_INTERPOLATION_PRESSURE,BASIS_XI_INTERPOLATION_PRESSURE, & 
466        & BASIS_XI_INTERPOLATION_PRESSURE/),Err)                         
467      IF(BASIS_TYPE/=CMISS_BASIS_SIMPLEX_TYPE) THEN
468        CALL CMISSBasis_QuadratureNumberOfGaussXiSet(BasisPressure,(/BASIS_GAUSS_PRESSURE,BASIS_GAUSS_PRESSURE, & 
469          & BASIS_GAUSS_PRESSURE/),Err)
470      ELSE
471        CALL CMISSBasis_QuadratureOrderSet(BasisPressure,BASIS_GAUSS_PRESSURE+1,Err)
472      ENDIF
473    ENDIF
474    !Finish the creation of the basis
475    CALL CMISSBasis_CreateFinish(BasisPressure,Err)
476  ENDIF
477
478  !
479  !================================================================================================================================
480  !
481
482  !MESH
483
484  !Start the creation of mesh nodes
485  CALL CMISSNodes_Initialise(Nodes,Err)
486  CALL CMISSNodes_CreateStart(Region,TOTAL_NUMBER_OF_NODES,Nodes,Err)
487  CALL CMISSNodes_CreateFinish(Nodes,Err)
488  !Start the creation of the mesh
489  CALL CMISSMesh_Initialise(Mesh,Err)
490  CALL CMISSMesh_CreateStart(MeshUserNumber,Region,NUMBER_OF_DIMENSIONS,Mesh,Err)
491  !Set number of mesh elements
492  CALL CMISSMesh_NumberOfElementsSet(Mesh,TOTAL_NUMBER_OF_ELEMENTS,Err)
493  !Set number of mesh components
494  CALL CMISSMesh_NumberOfComponentsSet(Mesh,MESH_NUMBER_OF_COMPONENTS,Err)
495  !Specify spatial mesh component
496  CALL CMISSMeshElements_Initialise(MeshElementsSpace,Err)
497  CALL CMISSMeshElements_Initialise(MeshElementsVelocity,Err)
498  CALL CMISSMeshElements_Initialise(MeshElementsPressure,Err)
499  MESH_COMPONENT_NUMBER_SPACE=1
500  MESH_COMPONENT_NUMBER_VELOCITY=1
501  MESH_COMPONENT_NUMBER_PRESSURE=1
502  CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_SPACE,BasisSpace,MeshElementsSpace,Err)
503  DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
504    CALL CMISSMeshElements_NodesSet(MeshElementsSpace,ELEMENT_NUMBER,CM%M(ELEMENT_NUMBER,1:NUMBER_OF_ELEMENT_NODES_SPACE),Err)
505  ENDDO
506  CALL CMISSMeshElements_CreateFinish(MeshElementsSpace,Err)
507  !Specify velocity mesh component
508  IF(BASIS_XI_INTERPOLATION_VELOCITY==BASIS_XI_INTERPOLATION_SPACE) THEN
509    MeshElementsVelocity=MeshElementsSpace
510  ELSE
511    MESH_COMPONENT_NUMBER_VELOCITY=MESH_COMPONENT_NUMBER_SPACE+1
512    CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_VELOCITY,BasisVelocity,MeshElementsVelocity,Err)
513    DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
514      CALL CMISSMeshElements_NodesSet(MeshElementsVelocity,ELEMENT_NUMBER,CM%V(ELEMENT_NUMBER, & 
515        & 1:NUMBER_OF_ELEMENT_NODES_VELOCITY),Err)
516    ENDDO
517    CALL CMISSMeshElements_CreateFinish(MeshElementsVelocity,Err)
518  ENDIF
519  !Specify pressure mesh component
520  IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_SPACE) THEN
521    MeshElementsPressure=MeshElementsSpace
522    MESH_COMPONENT_NUMBER_PRESSURE=MESH_COMPONENT_NUMBER_SPACE
523  ELSE IF(BASIS_XI_INTERPOLATION_PRESSURE==BASIS_XI_INTERPOLATION_VELOCITY) THEN
524    MeshElementsPressure=MeshElementsVelocity
525    MESH_COMPONENT_NUMBER_PRESSURE=MESH_COMPONENT_NUMBER_VELOCITY
526  ELSE
527    MESH_COMPONENT_NUMBER_PRESSURE=MESH_COMPONENT_NUMBER_VELOCITY+1
528    CALL CMISSMeshElements_CreateStart(Mesh,MESH_COMPONENT_NUMBER_PRESSURE,BasisPressure,MeshElementsPressure,Err)
529    DO ELEMENT_NUMBER=1,TOTAL_NUMBER_OF_ELEMENTS
530      CALL CMISSMeshElements_NodesSet(MeshElementsPressure,ELEMENT_NUMBER,CM%P(ELEMENT_NUMBER, & 
531        & 1:NUMBER_OF_ELEMENT_NODES_PRESSURE),Err)
532    ENDDO
533    CALL CMISSMeshElements_CreateFinish(MeshElementsPressure,Err)
534  ENDIF
535  !Finish the creation of the mesh
536  CALL CMISSMesh_CreateFinish(Mesh,Err)
537
538  !
539  !================================================================================================================================
540  !
541
542  !GEOMETRIC FIELD
543
544  !Create a decomposition
545  CALL CMISSDecomposition_Initialise(Decomposition,Err)
546  CALL CMISSDecomposition_CreateStart(DecompositionUserNumber,Mesh,Decomposition,Err)
547  !Set the decomposition to be a general decomposition with the specified number of domains
548  CALL CMISSDecomposition_TypeSet(Decomposition,CMISS_DECOMPOSITION_CALCULATED_TYPE,Err)
549  CALL CMISSDecomposition_NumberOfDomainsSet(Decomposition,DomainUserNumber,Err)
550  !Finish the decomposition
551  CALL CMISSDecomposition_CreateFinish(Decomposition,Err)
552
553  !Start to create a default (geometric) field on the region
554  CALL CMISSField_Initialise(GeometricField,Err)
555  CALL CMISSField_CreateStart(GeometricFieldUserNumber,Region,GeometricField,Err)
556  !Set the field type
557  CALL CMISSField_TypeSet(GeometricField,CMISS_FIELD_GEOMETRIC_TYPE,Err)
558  !Set the decomposition to use
559  CALL CMISSField_MeshDecompositionSet(GeometricField,Decomposition,Err)
560  !Set the scaling to use
561  CALL CMISSField_ScalingTypeSet(GeometricField,CMISS_FIELD_NO_SCALING,Err)
562  !Set the mesh component to be used by the field components.
563  DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
564    CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, & 
565      & MESH_COMPONENT_NUMBER_SPACE,Err)
566  ENDDO
567  !Finish creating the field
568  CALL CMISSField_CreateFinish(GeometricField,Err)
569  !Update the geometric field parameters
570  DO NODE_NUMBER=1,NUMBER_OF_NODES_SPACE
571    DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
572      VALUE=CM%N(NODE_NUMBER,COMPONENT_NUMBER)*0.3141592653589793_CMISSDP
573      CALL CMISSField_ParameterSetUpdateNode(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 
574        & 1,CMISS_NO_GLOBAL_DERIV,NODE_NUMBER,COMPONENT_NUMBER,VALUE,Err)
575    ENDDO
576  ENDDO
577  CALL CMISSField_ParameterSetUpdateStart(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,Err)
578  CALL CMISSField_ParameterSetUpdateFinish(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE,Err)
579
580  !
581  !================================================================================================================================
582  !
583
584  !EQUATIONS SETS
585
586  !Create the equations set for dynamic Navier-Stokes
587  CALL CMISSEquationsSet_Initialise(EquationsSetNavierStokes,Err)
588  CALL CMISSField_Initialise(EquationsSetField,Err)
589  CALL CMISSEquationsSet_CreateStart(EquationsSetUserNumberNavierStokes,Region,GeometricField, &
590    & CMISS_EQUATIONS_SET_FLUID_MECHANICS_CLASS,CMISS_EQUATIONS_SET_NAVIER_STOKES_EQUATION_TYPE,&
591    & CMISS_EQUATIONS_SET_TRANSIENT_NAVIER_STOKES_SUBTYPE, &
592    & EquationsSetFieldUserNumber,EquationsSetField,EquationsSetNavierStokes,Err)
593  !Set the equations set to be a dynamic Navier-Stokes problem
594  
595  !Finish creating the equations set
596  CALL CMISSEquationsSet_CreateFinish(EquationsSetNavierStokes,Err)
597
598
599  !
600  !================================================================================================================================
601  !
602
603  !DEPENDENT FIELDS
604
605  !Create the equations set dependent field variables for dynamic Navier-Stokes
606  CALL CMISSField_Initialise(DependentFieldNavierStokes,Err)
607  CALL CMISSEquationsSet_DependentCreateStart(EquationsSetNavierStokes,DependentFieldUserNumberNavierStokes, & 
608    & DependentFieldNavierStokes,Err)
609  !Set the mesh component to be used by the field components.
610  DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
611    CALL CMISSField_ComponentMeshComponentSet(DependentFieldNavierStokes,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, & 
612      & MESH_COMPONENT_NUMBER_VELOCITY,Err)
613    CALL CMISSField_ComponentMeshComponentSet(DependentFieldNavierStokes,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, & 
614      & MESH_COMPONENT_NUMBER_VELOCITY,Err)
615  ENDDO
616  COMPONENT_NUMBER=NUMBER_OF_DIMENSIONS+1
617    CALL CMISSField_ComponentMeshComponentSet(DependentFieldNavierStokes,CMISS_FIELD_U_VARIABLE_TYPE,COMPONENT_NUMBER, & 
618      & MESH_COMPONENT_NUMBER_PRESSURE,Err)
619    CALL CMISSField_ComponentMeshComponentSet(DependentFieldNavierStokes,CMISS_FIELD_DELUDELN_VARIABLE_TYPE,COMPONENT_NUMBER, & 
620      & MESH_COMPONENT_NUMBER_PRESSURE,Err)
621  !Finish the equations set dependent field variables
622  CALL CMISSEquationsSet_DependentCreateFinish(EquationsSetNavierStokes,Err)
623
624  !Initialise dependent field
625  DO COMPONENT_NUMBER=1,NUMBER_OF_DIMENSIONS
626    CALL CMISSField_ComponentValuesInitialise(DependentFieldNavierStokes,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 
627      & COMPONENT_NUMBER,INITIAL_FIELD_NAVIER_STOKES(COMPONENT_NUMBER),Err)
628  ENDDO
629
630
631  !
632  !================================================================================================================================
633  !
634
635  !MATERIALS FIELDS
636
637  !Create the equations set materials field variables for dynamic Navier-Stokes
638  CALL CMISSField_Initialise(MaterialsFieldNavierStokes,Err)
639  CALL CMISSEquationsSet_MaterialsCreateStart(EquationsSetNavierStokes,MaterialsFieldUserNumberNavierStokes, & 
640    & MaterialsFieldNavierStokes,Err)
641  !Finish the equations set materials field variables
642  CALL CMISSEquationsSet_MaterialsCreateFinish(EquationsSetNavierStokes,Err)
643  CALL CMISSField_ComponentValuesInitialise(MaterialsFieldNavierStokes,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 
644    & MaterialsFieldUserNumberNavierStokesMu,MU_PARAM_NAVIER_STOKES,Err)
645  CALL CMISSField_ComponentValuesInitialise(MaterialsFieldNavierStokes,CMISS_FIELD_U_VARIABLE_TYPE,CMISS_FIELD_VALUES_SET_TYPE, & 
646    & MaterialsFieldUserNumberNavierStokesRho,RHO_PARAM_NAVIER_STOKES,Err)
647
648  !
649  !================================================================================================================================
650  !
651
652  !ANALYTIC FIELDS
653
654  !Create the equations set analytic field variables for static Navier-Stokes
655  CALL CMISSField_Initialise(AnalyticFieldNavierStokes,Err)
656  IF(NUMBER_OF_DIMENSIONS==2) THEN  
657    CALL CMISSEquationsSet_AnalyticCreateStart(EquationsSetNavierStokes,ANALYTICAL_TYPE,AnalyticFieldUserNumberNavierStokes, &
658      & AnalyticFieldNavierStokes,Err)
659  ELSE
660    CALL CMISSEquationsSet_AnalyticCreateStart(EquationsSetNavierStokes,ANALYTICAL_TYPE,AnalyticFieldUserNumberNavierStokes, &
661      & AnalyticFieldNavierStokes,Err)
662  ENDIF
663  !Finish the equations set analytic field variables
664  CALL CMISSEquationsSet_AnalyticCreateFinish(EquationsSetNavierStokes,Err)
665
666  !
667  !================================================================================================================================
668  !
669
670  !EQUATIONS
671
672
673  !Create the equations set equations
674  CALL CMISSEquations_Initialise(EquationsNavierStokes,Err)
675  CALL CMISSEquationsSet_EquationsCreateStart(EquationsSetNavierStokes,EquationsNavierStokes,Err)
676  !Set the equations matrices sparsity type
677  CALL CMISSEquations_SparsityTypeSet(EquationsNavierStokes,CMISS_EQUATIONS_SPARSE_MATRICES,Err)
678  !Set the equations lumping type
679  CALL CMISSEquations_LumpingTypeSet(EquationsNavierStokes,CMISS_EQUATIONS_UNLUMPED_MATRICES,Err)
680  !Set the equations set output
681  CALL CMISSEquations_OutputTypeSet(EquationsNavierStokes,EQUATIONS_NAVIER_STOKES_OUTPUT,Err)
682  !Finish the equations set equations
683  CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSetNavierStokes,Err)
684
685
686  !
687  !================================================================================================================================
688  !
689
690  !PROBLEMS
691
692  !Start the creation of a problem.
693  CALL CMISSProblem_Initialise(Problem,Err)
694  CALL CMISSControlLoop_Initialise(ControlLoop,Err)
695  CALL CMISSProblem_CreateStart(ProblemUserNumber,Problem,Err)
696  !Set the problem to be a dynamic Navier-Stokes problem
697  CALL CMISSProblem_SpecificationSet(Problem,CMISS_PROBLEM_FLUID_MECHANICS_CLASS,CMISS_PROBLEM_NAVIER_STOKES_EQUATION_TYPE, &
698    & CMISS_PROBLEM_TRANSIENT_NAVIER_STOKES_SUBTYPE,Err)
699  !Finish the creation of a problem.
700  CALL CMISSProblem_CreateFinish(Problem,Err)
701  !Start the creation of the problem control loop
702  CALL CMISSProblem_ControlLoopCreateStart(Problem,Err)
703  !Get the control loop
704  CALL CMISSProblem_ControlLoopGet(Problem,CMISS_CONTROL_LOOP_NODE,ControlLoop,Err)
705  !Set the times
706  CALL CMISSControlLoop_TimesSet(ControlLoop,DYNAMIC_SOLVER_NAVIER_STOKES_START_TIME,DYNAMIC_SOLVER_NAVIER_STOKES_STOP_TIME, & 
707    & DYNAMIC_SOLVER_NAVIER_STOKES_TIME_INCREMENT,Err)
708  !Set the output timing
709  CALL CMISSControlLoop_TimeOutputSet(ControlLoop,DYNAMIC_SOLVER_NAVIER_STOKES_OUTPUT_FREQUENCY,Err)
710  !Finish creating the problem control loop
711  CALL CMISSProblem_ControlLoopCreateFinish(Problem,Err)
712
713  !
714  !================================================================================================================================
715  !
716
717  !SOLVERS
718
719  !Start the creation of the problem solvers
720  CALL CMISSSolver_Initialise(DynamicSolverNavierStokes,Err)
721  CALL CMISSSolver_Initialise(NonlinearSolverNavierStokes,Err)
722  CALL CMISSSolver_Initialise(LinearSolverNavierStokes,Err)
723  CALL CMISSProblem_SolversCreateStart(Problem,Err)
724  !Get the dynamic dymamic solver
725  CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,SolverNavierStokesUserNumber,DynamicSolverNavierStokes,Err)
726  !Set the output type
727  CALL CMISSSolver_OutputTypeSet(DynamicSolverNavierStokes,DYNAMIC_SOLVER_NAVIER_STOKES_OUTPUT_TYPE,Err)
728  !Set theta
729  CALL CMISSSolver_DynamicThetaSet(DynamicSolverNavierStokes,DYNAMIC_SOLVER_NAVIER_STOKES_THETA,Err)
730!   CALL CMISSSolverDynamicDynamicSet(DynamicSolverNavierStokes,.TRUE.,Err)
731  !Get the dynamic nonlinear solver
732  CALL CMISSSolver_DynamicNonlinearSolverGet(DynamicSolverNavierStokes,NonlinearSolverNavierStokes,Err)
733  !Set the nonlinear Jacobian type
734  CALL CMISSSolver_NewtonJacobianCalculationTypeSet(NonlinearSolverNavierStokes,CMISS_SOLVER_NEWTON_JACOBIAN_EQUATIONS_CALCULATED, &
735    & Err)
736  !Set the output type
737  CALL CMISSSolver_OutputTypeSet(NonlinearSolverNavierStokes,NONLINEAR_SOLVER_NAVIER_STOKES_OUTPUT_TYPE,Err)
738  !Set the solver settings
739  CALL CMISSSolver_NewtonAbsoluteToleranceSet(NonlinearSolverNavierStokes,ABSOLUTE_TOLERANCE,Err)
740  CALL CMISSSolver_NewtonRelativeToleranceSet(NonlinearSolverNavierStokes,ABSOLUTE_TOLERANCE,Err)
741  CALL CMISSSolver_NewtonSolutionToleranceSet(NonlinearSolverNavierStokes,ABSOLUTE_TOLERANCE,Err)
742  !Get the dynamic nonlinear linear solver
743  CALL CMISSSolver_NewtonLinearSolverGet(NonlinearSolverNavierStokes,LinearSolverNavierStokes,Err)
744  !Set the output type
745  CALL CMISSSolver_OutputTypeSet(LinearSolverNavierStokes,LINEAR_SOLVER_NAVIER_STOKES_OUTPUT_TYPE,Err)
746  !Set the solver settings
747  IF(LINEAR_SOLVER_NAVIER_STOKES_DIRECT_FLAG) THEN
748    CALL CMISSSolver_LinearTypeSet(LinearSolverNavierStokes,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err)
749    CALL CMISSSolver_LibraryTypeSet(LinearSolverNavierStokes,CMISS_SOLVER_MUMPS_LIBRARY,Err)
750  ELSE
751    CALL CMISSSolver_LinearTypeSet(LinearSolverNavierStokes,CMISS_SOLVER_LINEAR_ITERATIVE_SOLVE_TYPE,Err)
752! ! !     CALL CMISSSolver_LinearIterativeMaximumIterationsSet(LinearSolverNavierStokes,MAXIMUM_ITERATIONS,Err)
753    CALL CMISSSolver_LinearIterativeDivergenceToleranceSet(LinearSolverNavierStokes,DIVERGENCE_TOLERANCE,Err)
754    CALL CMISSSolver_LinearIterativeRelativeToleranceSet(LinearSolverNavierStokes,RELATIVE_TOLERANCE,Err)
755    CALL CMISSSolver_LinearIterativeAbsoluteToleranceSet(LinearSolverNavierStokes,ABSOLUTE_TOLERANCE,Err)
756    CALL CMISSSolver_LinearIterativeGMRESRestartSet(LinearSolverNavierStokes,RESTART_VALUE,Err)
757  ENDIF
758  !Finish the creation of the problem solver
759  CALL CMISSProblem_SolversCreateFinish(Problem,Err)
760
761  !
762  !================================================================================================================================
763  !
764
765  !SOLVER EQUATIONS
766
767  !Start the creation of the problem solver equations
768  CALL CMISSSolver_Initialise(DynamicSolverNavierStokes,Err)
769  CALL CMISSSolverEquations_Initialise(SolverEquationsNavierStokes,Err)
770  CALL CMISSProblem_SolverEquationsCreateStart(Problem,Err)
771  !Get the dynamic solver equations
772  CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,SolverNavierStokesUserNumber,DynamicSolverNavierStokes,Err)
773  CALL CMISSSolver_SolverEquationsGet(DynamicSolverNavierStokes,SolverEquationsNavierStokes,Err)
774  !Set the solver equations sparsity
775  CALL CMISSSolverEquations_SparsityTypeSet(SolverEquationsNavierStokes,CMISS_SOLVER_SPARSE_MATRICES,Err)
776  !Add in the equations set
777  CALL CMISSSolverEquations_EquationsSetAdd(SolverEquationsNavierStokes,EquationsSetNavierStokes,EquationsSetIndex,Err)
778  !Finish the creation of the problem solver equations
779  CALL CMISSProblem_SolverEquationsCreateFinish(Problem,Err)
780
781  !
782  !================================================================================================================================
783  !
784
785  !BOUNDARY CONDITIONS
786
787  !Set up the boundary conditions as per the analytic solution
788  CALL CMISSBoundaryConditions_Initialise(BoundaryConditionsNavierStokes,Err)
789  CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquationsNavierStokes,BoundaryConditionsNavierStokes,Err)
790  CALL CMISSSolverEquations_BoundaryConditionsAnalytic(SolverEquationsNavierStokes,Err)
791  CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquationsNavierStokes,Err)
792
793  !
794  !================================================================================================================================
795  !
796
797  !RUN SOLVERS
798
799  !Output Analytic analysis
800  CALL CMISSAnalyticAnalysisOutput(DependentFieldNavierStokes,OUTPUT_STRING,Err)
801
802
803  !Turn of PETSc error handling
804  !CALL PETSC_ERRORHANDLING_SET_ON(ERR,ERROR,*999)
805
806  !Solve the problem
807  WRITE(*,'(A)') "Solving problem..."
808  CALL CMISSProblem_Solve(Problem,Err)
809  WRITE(*,'(A)') "Problem solved!"
810
811  !
812  !================================================================================================================================
813  !
814
815  !OUTPUT
816
817!   !Output Analytic analysis
818!   CALL CMISSAnalyticAnalysisOutput(DependentFieldNavierStokes,OUTPUT_STRING,Err)
819
820
821  EXPORT_FIELD_IO=.TRUE.
822  IF(EXPORT_FIELD_IO) THEN
823    WRITE(*,'(A)') "Exporting fields..."
824    CALL CMISSFields_Initialise(Fields,Err)
825    CALL CMISSFields_Create(Region,Fields,Err)
826    CALL CMISSFields_NodesExport(Fields,"DynamicNavierStokes","FORTRAN",Err)
827    CALL CMISSFields_ElementsExport(Fields,"DynamicNavierStokes","FORTRAN",Err)
828    CALL CMISSFields_Finalise(Fields,Err)
829    WRITE(*,'(A)') "Field exported!"
830  ENDIF
831  
832  !Finialise CMISS
833!   CALL CMISSFinalise(Err)
834
835  WRITE(*,'(A)') "Program successfully completed."
836  
837  STOP
838
839END PROGRAM ANALYTICNAVIERSTOKESEXAMPLE