/ClassicalField/Laplace/LaplaceEllipsoid/src/LaplaceEllipsoidExample.f90
FORTRAN Modern | 360 lines | 190 code | 52 blank | 118 comment | 0 complexity | 98da8b8a503479908fe314a7dbe4de75 MD5 | raw file
1!> \file 2!> \author Chris Bradley 3!> \brief This is an example program to solve a Laplace 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 ClassicalField/Laplace/Laplace/src/NewLaplaceExample.f90 43!! Example program to solve a Laplace equation using OpenCMISS calls. 44!! \htmlinclude ClassicalField/Laplace/Laplace/history.html 45!! 46!< 47 48!> Main program 49PROGRAM LAPLACEELLIPSOIDEXAMPLE 50 51 USE OPENCMISS 52 USE MPI 53 54 55#ifdef WIN32 56 USE IFQWIN 57#endif 58 59 IMPLICIT NONE 60 61 INTEGER(CMISSIntg), PARAMETER :: EquationsSetFieldUserNumber=1337 62 TYPE(CMISSFieldType) :: EquationsSetField 63 64 65 !Test program parameters 66 67 REAL(CMISSDP), PARAMETER :: LONG_AXIS=2.0_CMISSDP 68 REAL(CMISSDP), PARAMETER :: SHORT_AXIS=1.0_CMISSDP 69 REAL(CMISSDP), PARAMETER :: WALL_THICKNESS=0.5_CMISSDP 70 REAL(CMISSDP), PARAMETER :: CUTOFF_ANGLE=1.5708_CMISSDP 71 72 INTEGER(CMISSIntg), PARAMETER :: CoordinateSystemUserNumber=1 73 INTEGER(CMISSIntg), PARAMETER :: RegionUserNumber=2 74 INTEGER(CMISSIntg), PARAMETER :: QuadraticBasisUserNumber=1 75 INTEGER(CMISSIntg), PARAMETER :: QuadraticCollapsedBasisUserNumber=2 76 INTEGER(CMISSIntg), PARAMETER :: GeneratedMeshUserNumber=4 77 INTEGER(CMISSIntg), PARAMETER :: MeshUserNumber=5 78 INTEGER(CMISSIntg), PARAMETER :: DecompositionUserNumber=6 79 INTEGER(CMISSIntg), PARAMETER :: GeometricFieldUserNumber=7 80 INTEGER(CMISSIntg), PARAMETER :: DependentFieldUserNumber=8 81 INTEGER(CMISSIntg), PARAMETER :: EquationsSetUserNumber=9 82 INTEGER(CMISSIntg), PARAMETER :: ProblemUserNumber=10 83 INTEGER(CMISSIntg), PARAMETER :: NumberOfXiCoordinates=3 84!Program types 85 86 !Program variables 87 88 INTEGER(CMISSIntg) :: NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS,NUMBER_GLOBAL_Z_ELEMENTS 89 INTEGER(CMISSIntg) :: NUMBER_OF_DOMAINS 90 91 INTEGER(CMISSIntg) :: MPI_IERROR 92 93 LOGICAL :: EXPORT_FIELD 94 95 !CMISS variables 96 97 TYPE(CMISSBasisType) :: QuadraticBasis,QuadraticCollapsedBasis 98 TYPE(CMISSBoundaryConditionsType) :: BoundaryConditions 99 TYPE(CMISSCoordinateSystemType) :: CoordinateSystem,WorldCoordinateSystem 100 TYPE(CMISSDecompositionType) :: Decomposition 101 TYPE(CMISSEquationsType) :: Equations 102 TYPE(CMISSEquationsSetType) :: EquationsSet 103 TYPE(CMISSFieldType) :: GeometricField,DependentField 104 TYPE(CMISSFieldsType) :: Fields 105 TYPE(CMISSGeneratedMeshType) :: GeneratedMesh 106 TYPE(CMISSMeshType) :: Mesh 107 TYPE(CMISSProblemType) :: Problem 108 TYPE(CMISSRegionType) :: Region,WorldRegion 109 TYPE(CMISSSolverType) :: Solver 110 TYPE(CMISSSolverEquationsType) :: SolverEquations 111 112#ifdef WIN32 113 !Quickwin type 114 LOGICAL :: QUICKWIN_STATUS=.FALSE. 115 TYPE(WINDOWCONFIG) :: QUICKWIN_WINDOW_CONFIG 116#endif 117 118 !Generic CMISS variables 119 120 INTEGER(CMISSIntg) :: NumberOfComputationalNodes,ComputationalNodeNumber 121 INTEGER(CMISSIntg) :: EquationsSetIndex 122 INTEGER(CMISSIntg) :: FirstNodeNumber,LastNodeNumber 123 INTEGER(CMISSIntg) :: FirstNodeDomain,LastNodeDomain 124 INTEGER(CMISSIntg) :: Err 125 126#ifdef WIN32 127 !Initialise QuickWin 128 QUICKWIN_WINDOW_CONFIG%TITLE="General Output" !Window title 129 QUICKWIN_WINDOW_CONFIG%NUMTEXTROWS=-1 !Max possible number of rows 130 QUICKWIN_WINDOW_CONFIG%MODE=QWIN$SCROLLDOWN 131 !Set the window parameters 132 QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG) 133 !If attempt fails set with system estimated values 134 IF(.NOT.QUICKWIN_STATUS) QUICKWIN_STATUS=SETWINDOWCONFIG(QUICKWIN_WINDOW_CONFIG) 135#endif 136 137 !Intialise OpenCMISS 138 CALL CMISSInitialise(WorldCoordinateSystem,WorldRegion,Err) 139 140 CALL CMISSErrorHandlingModeSet(CMISS_ERRORS_TRAP_ERROR,Err) 141 142 !CALL CMISSDiagnosticsSetOn(CMISS_ALL_DIAG_TYPE,(/1,2,3,4,5/),"Diagnostics",(/"FIELD_MAPPINGS_CALCULATE"/),Err) 143 144 !Get the computational nodes information 145 CALL CMISSComputationalNumberOfNodesGet(NumberOfComputationalNodes,Err) 146 CALL CMISSComputationalNodeNumberGet(ComputationalNodeNumber,Err) 147 148 NUMBER_GLOBAL_X_ELEMENTS=4 149 NUMBER_GLOBAL_Y_ELEMENTS=3 150 NUMBER_GLOBAL_Z_ELEMENTS=1 151 NUMBER_OF_DOMAINS=NumberOfComputationalNodes 152 153 !Broadcast the number of elements in the X & Y directions and the number of partitions to the other computational nodes 154 CALL MPI_BCAST(NUMBER_GLOBAL_X_ELEMENTS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR) 155 CALL MPI_BCAST(NUMBER_GLOBAL_Y_ELEMENTS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR) 156 CALL MPI_BCAST(NUMBER_GLOBAL_Z_ELEMENTS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR) 157 CALL MPI_BCAST(NUMBER_OF_DOMAINS,1,MPI_INTEGER,0,MPI_COMM_WORLD,MPI_IERROR) 158 159 !Start the creation of a new RC coordinate system 160 CALL CMISSCoordinateSystem_Initialise(CoordinateSystem,Err) 161 CALL CMISSCoordinateSystem_CreateStart(CoordinateSystemUserNumber,CoordinateSystem,Err) 162 IF(NUMBER_GLOBAL_Z_ELEMENTS==0) THEN 163 !Set the coordinate system to be 2D 164 CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,2,Err) 165 ELSE 166 !Set the coordinate system to be 3D 167 CALL CMISSCoordinateSystem_DimensionSet(CoordinateSystem,3,Err) 168 ENDIF 169 !Finish the creation of the coordinate system 170 CALL CMISSCoordinateSystem_CreateFinish(CoordinateSystem,Err) 171 172 !Start the creation of the region 173 CALL CMISSRegion_Initialise(Region,Err) 174 CALL CMISSRegion_CreateStart(RegionUserNumber,WorldRegion,Region,Err) 175 !Set the regions coordinate system to the 2D RC coordinate system that we have created 176 CALL CMISSRegion_CoordinateSystemSet(Region,CoordinateSystem,Err) 177 !Finish the creation of the region 178 CALL CMISSRegion_CreateFinish(Region,Err) 179 180 !Define basis functions - tri-linear Lagrange and tri-Quadratic Lagrange, each with collapsed variant 181 !Quadratic Basis 182 CALL CMISSBasis_Initialise(QuadraticBasis,Err) 183 CALL CMISSBasis_CreateStart(QuadraticBasisUserNumber,QuadraticBasis,Err) 184 CALL CMISSBasis_InterpolationXiSet(QuadraticBasis,(/CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION, & 185 & CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION,CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION/),Err) 186 CALL CMISSBasis_QuadratureNumberOfGaussXiSet(QuadraticBasis, & 187 & (/CMISS_BASIS_MID_QUADRATURE_SCHEME,CMISS_BASIS_MID_QUADRATURE_SCHEME,CMISS_BASIS_MID_QUADRATURE_SCHEME/),Err) 188 CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(QuadraticBasis,.true.,Err) !Have to do this 189 CALL CMISSBasis_CreateFinish(QuadraticBasis,Err) 190 191 !Collapsed Quadratic Basis 192 CALL CMISSBasis_Initialise(QuadraticCollapsedBasis,Err) 193 CALL CMISSBasis_CreateStart(QuadraticCollapsedBasisUserNumber,QuadraticCollapsedBasis,Err) 194 CALL CMISSBasis_TypeSet(QuadraticCollapsedBasis,CMISS_BASIS_LAGRANGE_HERMITE_TP_TYPE,Err) 195 CALL CMISSBasis_NumberOfXiSet(QuadraticCollapsedBasis,NumberOfXiCoordinates,Err) 196 CALL CMISSBasis_InterpolationXiSet(QuadraticCollapsedBasis,(/CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION, & 197 & CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION,CMISS_BASIS_QUADRATIC_LAGRANGE_INTERPOLATION/),Err) 198 CALL CMISSBasis_CollapsedXiSet(QuadraticCollapsedBasis,(/CMISS_BASIS_XI_COLLAPSED, & 199 & CMISS_BASIS_COLLAPSED_AT_XI0,CMISS_BASIS_NOT_COLLAPSED/),Err) 200 CALL CMISSBasis_QuadratureNumberOfGaussXiSet(QuadraticCollapsedBasis, & 201 & (/CMISS_BASIS_MID_QUADRATURE_SCHEME,CMISS_BASIS_MID_QUADRATURE_SCHEME,CMISS_BASIS_MID_QUADRATURE_SCHEME/),Err) 202 CALL CMISSBasis_QuadratureLocalFaceGaussEvaluateSet(QuadraticCollapsedBasis,.true.,Err) !Have to do this 203 CALL CMISSBasis_CreateFinish(QuadraticCollapsedBasis,Err) 204 205 !Start the creation of a generated mesh in the region 206 CALL CMISSGeneratedMesh_Initialise(GeneratedMesh,Err) 207 CALL CMISSGeneratedMesh_CreateStart(GeneratedMeshUserNumber,Region,GeneratedMesh,Err) 208 !Set up a regular x*y*z mesh 209 CALL CMISSGeneratedMesh_TypeSet(GeneratedMesh,CMISS_GENERATED_MESH_ELLIPSOID_MESH_TYPE,Err) 210 !Set the default basis 211 CALL CMISSGeneratedMesh_BasisSet(GeneratedMesh,[QuadraticBasis,QuadraticCollapsedBasis],Err) 212 !Define the mesh on the region 213 214 CALL CMISSGeneratedMesh_ExtentSet(GeneratedMesh,(/LONG_AXIS,SHORT_AXIS,WALL_THICKNESS,CUTOFF_ANGLE/),Err) 215 CALL CMISSGeneratedMesh_NumberOfElementsSet(GeneratedMesh,(/NUMBER_GLOBAL_X_ELEMENTS,NUMBER_GLOBAL_Y_ELEMENTS, & 216 & NUMBER_GLOBAL_Z_ELEMENTS/),Err) 217 218 !Finish the creation of a generated mesh in the region 219 CALL CMISSMesh_Initialise(Mesh,Err) 220 CALL CMISSGeneratedMesh_CreateFinish(GeneratedMesh,MeshUserNumber,Mesh,Err) 221 222 !Create a decomposition 223 CALL CMISSDecomposition_Initialise(Decomposition,Err) 224 CALL CMISSDecomposition_CreateStart(DecompositionUserNumber,Mesh,Decomposition,Err) 225 !Set the decomposition to be a general decomposition with the specified number of domains 226 CALL CMISSDecomposition_TypeSet(Decomposition,CMISS_DECOMPOSITION_CALCULATED_TYPE,Err) 227 CALL CMISSDecomposition_NumberOfDomainsSet(Decomposition,NUMBER_OF_DOMAINS,Err) 228 !Finish the decomposition 229 CALL CMISSDecomposition_CreateFinish(Decomposition,Err) 230 231 !Start to create a default (geometric) field on the region 232 CALL CMISSField_Initialise(GeometricField,Err) 233 CALL CMISSField_CreateStart(GeometricFieldUserNumber,Region,GeometricField,Err) 234 !Set the decomposition to use 235 CALL CMISSField_MeshDecompositionSet(GeometricField,Decomposition,Err) 236 !Set the domain to be used by the field components. 237 CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,1,1,Err) 238 CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,2,1,Err) 239 IF(NUMBER_GLOBAL_Z_ELEMENTS/=0) THEN 240 CALL CMISSField_ComponentMeshComponentSet(GeometricField,CMISS_FIELD_U_VARIABLE_TYPE,3,1,Err) 241 ENDIF 242 !Finish creating the field 243 CALL CMISSField_CreateFinish(GeometricField,Err) 244 245 !Update the geometric field parameters 246 CALL CMISSGeneratedMesh_GeometricParametersCalculate(GeneratedMesh,GeometricField,Err) 247 248 !Create the equations_set 249 CALL CMISSEquationsSet_Initialise(EquationsSet,Err) 250 CALL CMISSField_Initialise(EquationsSetField,Err) 251CALL CMISSEquationsSet_CreateStart(EquationsSetUserNumber,Region,GeometricField,CMISS_EQUATIONS_SET_CLASSICAL_FIELD_CLASS, & 252 & CMISS_EQUATIONS_SET_LAPLACE_EQUATION_TYPE,CMISS_EQUATIONS_SET_STANDARD_LAPLACE_SUBTYPE,EquationsSetFieldUserNumber, & 253 & EquationsSetField, & 254 & EquationsSet,Err) 255 !Set the equations set to be a standard Laplace problem 256 257 !Finish creating the equations set 258 CALL CMISSEquationsSet_CreateFinish(EquationsSet,Err) 259 260 !Create the equations set dependent field variables 261 CALL CMISSField_Initialise(DependentField,Err) 262 CALL CMISSEquationsSet_DependentCreateStart(EquationsSet,DependentFieldUserNumber,DependentField,Err) 263 !Finish the equations set dependent field variables 264 CALL CMISSEquationsSet_DependentCreateFinish(EquationsSet,Err) 265 266 !Create the equations set equations 267 CALL CMISSEquations_Initialise(Equations,Err) 268 CALL CMISSEquationsSet_EquationsCreateStart(EquationsSet,Equations,Err) 269 !Set the equations matrices sparsity type 270 CALL CMISSEquations_SparsityTypeSet(Equations,CMISS_EQUATIONS_SPARSE_MATRICES,Err) 271 !Set the equations set output 272 !CALL CMISSEquations_OutputTypeSet(Equations,CMISS_EQUATIONS_NO_OUTPUT,Err) 273 CALL CMISSEquations_OutputTypeSet(Equations,CMISS_EQUATIONS_TIMING_OUTPUT,Err) 274 !CALL CMISSEquations_OutputTypeSet(Equations,CMISS_EQUATIONS_MATRIX_OUTPUT,Err) 275 !CALL CMISSEquations_OutputTypeSet(Equations,CMISS_EQUATIONS_ELEMENT_MATRIX_OUTPUT,Err) 276 !Finish the equations set equations 277 CALL CMISSEquationsSet_EquationsCreateFinish(EquationsSet,Err) 278 279 !Start the creation of a problem. 280 CALL CMISSProblem_Initialise(Problem,Err) 281 CALL CMISSProblem_CreateStart(ProblemUserNumber,Problem,Err) 282 !Set the problem to be a standard Laplace problem 283 CALL CMISSProblem_SpecificationSet(Problem,CMISS_PROBLEM_CLASSICAL_FIELD_CLASS,CMISS_PROBLEM_LAPLACE_EQUATION_TYPE, & 284 & CMISS_PROBLEM_STANDARD_LAPLACE_SUBTYPE,Err) 285 !Finish the creation of a problem. 286 CALL CMISSProblem_CreateFinish(Problem,Err) 287 288 !Start the creation of the problem control loop 289 CALL CMISSProblem_ControlLoopCreateStart(Problem,Err) 290 !Finish creating the problem control loop 291 CALL CMISSProblem_ControlLoopCreateFinish(Problem,Err) 292 293 !Start the creation of the problem solvers 294 CALL CMISSSolver_Initialise(Solver,Err) 295 CALL CMISSProblem_SolversCreateStart(Problem,Err) 296 CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,1,Solver,Err) 297 !CALL CMISSSolver_OutputTypeSet(Solver,CMISS_SOLVER_NO_OUTPUT,Err) 298 !CALL CMISSSolver_OutputTypeSet(Solver,CMISS_SOLVER_PROGRESS_OUTPUT,Err) 299 !CALL CMISSSolver_OutputTypeSet(Solver,CMISS_SOLVER_TIMING_OUTPUT,Err) 300 !CALL CMISSSolver_OutputTypeSet(Solver,CMISS_SOLVER_SOLVER_OUTPUT,Err) 301 CALL CMISSSolver_OutputTypeSet(Solver,CMISS_SOLVER_MATRIX_OUTPUT,Err) 302 CALL CMISSSolver_LinearTypeSet(Solver,CMISS_SOLVER_LINEAR_DIRECT_SOLVE_TYPE,Err) 303 CALL CMISSSolver_LibraryTypeSet(Solver,CMISS_SOLVER_MUMPS_LIBRARY,Err) 304 !Finish the creation of the problem solver 305 CALL CMISSProblem_SolversCreateFinish(Problem,Err) 306 307 !Start the creation of the problem solver equations 308 CALL CMISSSolver_Initialise(Solver,Err) 309 CALL CMISSSolverEquations_Initialise(SolverEquations,Err) 310 CALL CMISSProblem_SolverEquationsCreateStart(Problem,Err) 311 !Get the solve equations 312 CALL CMISSProblem_SolverGet(Problem,CMISS_CONTROL_LOOP_NODE,1,Solver,Err) 313 CALL CMISSSolver_SolverEquationsGet(Solver,SolverEquations,Err) 314 !Set the solver equations sparsity 315 CALL CMISSSolverEquations_SparsityTypeSet(SolverEquations,CMISS_SOLVER_SPARSE_MATRICES,Err) 316 !CALL CMISSSolverEquations_SparsityTypeSet(SolverEquations,CMISS_SOLVER_FULL_MATRICES,Err) 317 !Add in the equations set 318 CALL CMISSSolverEquations_EquationsSetAdd(SolverEquations,EquationsSet,EquationsSetIndex,Err) 319 !Finish the creation of the problem solver equations 320 CALL CMISSProblem_SolverEquationsCreateFinish(Problem,Err) 321 322 !Start the creation of the equations set boundary conditions 323 CALL CMISSBoundaryConditions_Initialise(BoundaryConditions,Err) 324 CALL CMISSSolverEquations_BoundaryConditionsCreateStart(SolverEquations,BoundaryConditions,Err) 325 !Set the first node to 0.0 and the last node to 1.0 326 FirstNodeNumber=1 327 LastNodeNumber=2!SIZE(Region%region%Nodes) 328 !CALL CMISSDecomposition_NodeDomainGet(Decomposition,FirstNodeNumber,1,FirstNodeDomain,Err) 329 !CALL CMISSDecomposition_NodeDomainGet(Decomposition,LastNodeNumber,1,LastNodeDomain,Err) 330 IF(FirstNodeDomain==ComputationalNodeNumber) THEN 331 CALL CMISSBoundaryConditions_SetNode(BoundaryConditions,DependentField,CMISS_FIELD_U_VARIABLE_TYPE,1,1,FirstNodeNumber,1, & 332 & CMISS_BOUNDARY_CONDITION_FIXED,0.0_CMISSDP,Err) 333 ENDIF 334 IF(LastNodeDomain==ComputationalNodeNumber) THEN 335 CALL CMISSBoundaryConditions_SetNode(BoundaryConditions,DependentField,CMISS_FIELD_U_VARIABLE_TYPE,1,1,LastNodeNumber,1, & 336 & CMISS_BOUNDARY_CONDITION_FIXED,1.0_CMISSDP,Err) 337 ENDIF 338 !Finish the creation of the equations set boundary conditions 339 CALL CMISSSolverEquations_BoundaryConditionsCreateFinish(SolverEquations,Err) 340 341 !Solve the problem 342 CALL CMISSProblem_Solve(Problem,Err) 343 344 EXPORT_FIELD=.TRUE. 345 IF(EXPORT_FIELD) THEN 346 CALL CMISSFields_Initialise(Fields,Err) 347 CALL CMISSFields_Create(Region,Fields,Err) 348 CALL CMISSFields_NodesExport(Fields,"NewLaplace","FORTRAN",Err) 349 CALL CMISSFields_ElementsExport(Fields,"NewLaplace","FORTRAN",Err) 350 CALL CMISSFields_Finalise(Fields,Err) 351 ENDIF 352 353 !Finialise CMISS 354 CALL CMISSFinalise(Err) 355 356 WRITE(*,'(A)') "Program successfully completed." 357 358 STOP 359 360END PROGRAM LAPLACEELLIPSOIDEXAMPLE