fluxes_kt.f90
Go to the documentation of this file.
1!#############################################################################
2!# #
3!# fosite - 3D hydrodynamical simulation program #
4!# module: fluxes_kt.f90 #
5!# #
6!# Copyright (C) 2007-2018 #
7!# Tobias Illenseer <tillense@astrophysik.uni-kiel.de> #
8!# Manuel Jung <mjung@astrophysik.uni-kiel.de> #
9!# Jannes Klee <jklee@astrophysik.uni-kiel.de> #
10!# Jubin Lirawi <jlirawi@astrophysik.uni-kiel.de> #
11!# #
12!# This program is free software; you can redistribute it and/or modify #
13!# it under the terms of the GNU General Public License as published by #
14!# the Free Software Foundation; either version 2 of the License, or (at #
15!# your option) any later version. #
16!# #
17!# This program is distributed in the hope that it will be useful, but #
18!# WITHOUT ANY WARRANTY; without even the implied warranty of #
19!# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, GOOD TITLE or #
20!# NON INFRINGEMENT. See the GNU General Public License for more #
21!# details. #
22!# #
23!# You should have received a copy of the GNU General Public License #
24!# along with this program; if not, write to the Free Software #
25!# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. #
26!# #
27!#############################################################################
28
29!----------------------------------------------------------------------------!
39!----------------------------------------------------------------------------!
44!#ifdef PARALLEL
45! DEFAULT_MPI_REAL
46!#endif
49 USE common_dict
50#ifdef PARALLEL
51#ifdef HAVE_MPI_MOD
52 USE mpi
53#endif
54#endif
55 IMPLICIT NONE
56#ifdef PARALLEL
57#ifdef HAVE_MPIF_H
58 include 'mpif.h'
59#endif
60#endif
61 !--------------------------------------------------------------------------!
62 PRIVATE
63 TYPE, EXTENDS (fluxes_base) :: fluxes_kt
64 CONTAINS
65 PROCEDURE :: initfluxes_kt
66 PROCEDURE :: calculatefluxes
67 PROCEDURE :: finalize
68 END TYPE
69 !--------------------------------------------------------------------------!
70 PUBLIC :: &
71 ! types
73 !--------------------------------------------------------------------------!
74
75CONTAINS
76
77 SUBROUTINE initfluxes_kt(this,Mesh,Physics,config,IO)
79 IMPLICIT NONE
80 !------------------------------------------------------------------------!
81 CLASS(fluxes_kt), INTENT(INOUT) :: this
82 CLASS(mesh_base), INTENT(IN) :: Mesh
83 CLASS(physics_base), INTENT(IN) :: Physics
84 TYPE(dict_typ), INTENT(IN), POINTER :: config,IO
85 !------------------------------------------------------------------------!
86 CALL this%InitFluxes(mesh,physics,config,io,kt,"KT")
87
89 SELECT TYPE(mesh)
90 TYPE IS(mesh_midpoint)
91 ! do nothing
92! CASE(TRAPEZOIDAL)
93! CALL this%Error("InitFluxes_kt", "Trapezoidal not supported from > 0.6")
94 CLASS DEFAULT
95 CALL this%Error("InitFluxes_kt", "Unknown mesh type.")
96 END SELECT
97 END SUBROUTINE initfluxes_kt
98
99
109 SUBROUTINE calculatefluxes(this,Mesh,Physics,pvar,cvar, &
110 xfluxdydz,yfluxdzdx,zfluxdxdy)
111 IMPLICIT NONE
112 !------------------------------------------------------------------------!
113 CLASS(fluxes_kt), INTENT(INOUT) :: this
114 CLASS(mesh_base), INTENT(IN) :: Mesh
115 CLASS(physics_base), INTENT(INOUT) :: Physics
116 CLASS(marray_compound), INTENT(INOUT) :: pvar,cvar
117 REAL, DIMENSION(Mesh%IGMIN:Mesh%IGMAX,Mesh%JGMIN:Mesh%JGMAX,Mesh%KGMIN:Mesh%KGMAX,Physics%vnum), &
118 INTENT(OUT) :: xfluxdydz,yfluxdzdx,zfluxdxdy
119 !------------------------------------------------------------------------!
120 INTEGER :: i,j,k,l,m
121 !------------------------------------------------------------------------!
122 ! execute generic tasks common to all flux types
123 CALL this%CalculateFaceData(mesh,physics,pvar,cvar)
124
125 ! compute numerical fluxes along x-direction (west and east) devided by dydz
126 IF (mesh%INUM.GT.1) THEN
127 m = 1 ! index counting the number of spatial directions for which transport is enabled
128 ! physical fluxes
129 CALL physics%CalculateFluxesX(mesh,2*m-1,2*m,this%prim,this%cons,this%pfluxes)
130!NEC$ SHORTLOOP
131 DO l=1,physics%VNUM
132!NEC$ IVDEP
133 DO k=mesh%KGMIN,mesh%KGMAX
134!NEC$ IVDEP
135 DO j=mesh%JGMIN,mesh%JGMAX
136!NEC$ IVDEP
137 DO i=mesh%IMIN+mesh%IM1,mesh%IMAX
138 xfluxdydz(i,j,k,l) = mesh%dAxdydz(i+1,j,k,1) / &
139 (this%maxwav%data4d(i,j,k,m) - this%minwav%data4d(i,j,k,m)) * &
140 (this%maxwav%data4d(i,j,k,m) * this%pfluxes%data5d(i,j,k,2*m,l) - &
141 this%minwav%data4d(i,j,k,m) * this%pfluxes%data5d(i+1,j,k,2*m-1,l) + &
142 this%minwav%data4d(i,j,k,m) * this%maxwav%data4d(i,j,k,m) * &
143 (this%cons%data5d(i+1,j,k,2*m-1,l) - this%cons%data5d(i,j,k,2*m,l)))
144 END DO
145 END DO
146 END DO
147 END DO
148 ELSE
149 m = 0
150 xfluxdydz(:,:,:,:) = 0.0
151 END IF
152
153 ! compute numerical fluxes along y-direction (south and north) devided by dzdx
154 IF (mesh%JNUM.GT.1) THEN
155 m = m + 1 ! increase wave speed index, may be 1 or 2 now depending
156 ! on whether there was transport in x-direction or not
157 ! physical fluxes
158 CALL physics%CalculateFluxesY(mesh,2*m-1,2*m,this%prim,this%cons,this%pfluxes)
159!NEC$ SHORTLOOP
160 DO l=1,physics%VNUM
161!NEC$ IVDEP
162 DO k=mesh%KGMIN,mesh%KGMAX
163!NEC$ IVDEP
164 DO j=mesh%JMIN+mesh%JM1,mesh%JMAX
165!NEC$ IVDEP
166 DO i=mesh%IGMIN,mesh%IGMAX
167 yfluxdzdx(i,j,k,l) = mesh%dAydzdx(i,j+1,k,1) / &
168 (this%maxwav%data4d(i,j,k,m) - this%minwav%data4d(i,j,k,m)) * &
169 (this%maxwav%data4d(i,j,k,m) * this%pfluxes%data5d(i,j,k,2*m,l) - &
170 this%minwav%data4d(i,j,k,m) * this%pfluxes%data5d(i,j+1,k,2*m-1,l) + &
171 this%minwav%data4d(i,j,k,m) * this%maxwav%data4d(i,j,k,m) * &
172 (this%cons%data5d(i,j+1,k,2*m-1,l) - this%cons%data5d(i,j,k,2*m,l)))
173 END DO
174 END DO
175 END DO
176 END DO
177 ELSE
178 yfluxdzdx(:,:,:,:) = 0.0
179 END IF
180
181 ! compute numerical fluxes along z-direction (bottom and top) devided by dxdy
182 IF (mesh%KNUM.GT.1) THEN
183 m = m + 1 ! increase wave speed index, may be 1, 2 or 3 now depending
184 ! on whether there was transport in x- and/or y-direction or not
185 ! physical fluxes
186 CALL physics%CalculateFluxesZ(mesh,2*m-1,2*m,this%prim,this%cons,this%pfluxes)
187!NEC$ SHORTLOOP
188 DO l=1,physics%VNUM
189!NEC$ IVDEP
190 DO k=mesh%KMIN+mesh%KM1,mesh%KMAX
191!NEC$ IVDEP
192 DO j=mesh%JGMIN,mesh%JGMAX
193!NEC$ IVDEP
194 DO i=mesh%IGMIN,mesh%IGMAX
195 zfluxdxdy(i,j,k,l) = mesh%dAzdxdy(i,j,k+1,1) / &
196 (this%maxwav%data4d(i,j,k,m) - this%minwav%data4d(i,j,k,m)) * &
197 (this%maxwav%data4d(i,j,k,m) * this%pfluxes%data5d(i,j,k,2*m,l) - &
198 this%minwav%data4d(i,j,k,m) * this%pfluxes%data5d(i,j,k+1,2*m-1,l) + &
199 this%minwav%data4d(i,j,k,m) * this%maxwav%data4d(i,j,k,m) * &
200 (this%cons%data5d(i,j,k+1,2*m-1,l) - this%cons%data5d(i,j,k,2*m,l)))
201 END DO
202 END DO
203 END DO
204 END DO
205 ELSE
206 zfluxdxdy(:,:,:,:) = 0.0
207 END IF
208 END SUBROUTINE calculatefluxes
209
210 ! Here was once a routine for BilinearInterpolation, which was needed for
211 ! trapezoidal rule. See (<= 0.6) version of fosite.
212
213 SUBROUTINE finalize(this)
214 IMPLICIT NONE
215 !------------------------------------------------------------------------!
216 CLASS(fluxes_kt), INTENT(INOUT) :: this
217 !------------------------------------------------------------------------!
218 CALL this%Finalize_base()
219 END SUBROUTINE finalize
220
221END MODULE fluxes_kt_mod
Dictionary for generic data types.
Definition: common_dict.f90:61
base module for numerical flux functions
Definition: fluxes_base.f90:39
integer, parameter, public kt
numerical fluxes module for midpoint quadrature rule
Definition: fluxes_kt.f90:40
subroutine initfluxes_kt(this, Mesh, Physics, config, IO)
Definition: fluxes_kt.f90:78
subroutine calculatefluxes(this, Mesh, Physics, pvar, cvar, xfluxdydz, yfluxdzdx, zfluxdxdy)
Calculates the fluxes with the midpoint method.
Definition: fluxes_kt.f90:111
subroutine finalize(this)
Definition: fluxes_kt.f90:214
derived class for compound of mesh arrays
basic mesh module
Definition: mesh_base.f90:72
mesh module for midpoint quadrature rule
Basic physics module.
base module for reconstruction process
mesh data structure
Definition: mesh_base.f90:122