63 CLASS(reconstruction_base),
ALLOCATABLE &
73 REAL,
DIMENSION(:,:,:,:),
POINTER,
CONTIGUOUS &
79 procedure(calculatefluxes),
DEFERRED :: calculatefluxes
87 SUBROUTINE calculatefluxes(this,Mesh,Physics,pvar,cvar, &
88 xfluxdydz,yfluxdzdx,zfluxdxdy)
95 REAL,
DIMENSION(Mesh%IGMIN:Mesh%IGMAX,Mesh%JGMIN:Mesh%JGMAX,Mesh%KGMIN:Mesh%KGMAX,Physics%VNUM), &
96 INTENT(OUT) :: xfluxdydz,yfluxdzdx,zfluxdxdy
104 INTEGER,
PARAMETER ::
kt = 1
119 SUBROUTINE initfluxes(this,Mesh,Physics,config,IO,ftype,fname)
128 CHARACTER(LEN=*) :: fname
130 TYPE(
dict_typ),
POINTER :: IOrec => null()
131 INTEGER :: valwrite,i
132 CHARACTER(LEN=60) :: key
134 CALL this%InitLogging(ftype,fname)
137 IF (.NOT.mesh%Initialized().OR..NOT.physics%Initialized()) &
138 CALL this%Error(
"InitFluxes",
"mesh and/or physics module uninitialized")
141 ALLOCATE(this%minwav,this%maxwav, &
142 this%bxflux(mesh%JGMIN:mesh%JGMAX,mesh%KGMIN:mesh%KGMAX,2,physics%VNUM), &
143 this%byflux(mesh%KGMIN:mesh%KGMAX,mesh%IGMIN:mesh%IGMAX,2,physics%VNUM), &
144 this%bzflux(mesh%IGMIN:mesh%IGMAX,mesh%JGMIN:mesh%JGMAX,2,physics%VNUM), &
145 this%bxfold(mesh%JGMIN:mesh%JGMAX,mesh%KGMIN:mesh%KGMAX,2,physics%VNUM), &
146 this%byfold(mesh%KGMIN:mesh%KGMAX,mesh%IGMIN:mesh%IGMAX,2,physics%VNUM), &
147 this%bzfold(mesh%IGMIN:mesh%IGMAX,mesh%JGMIN:mesh%JGMAX,2,physics%VNUM), &
150 CALL this%Error(
"InitFluxes",
"Unable to allocate memory.")
158 CALL physics%new_statevector(this%prim,
primitive,mesh%NFACES)
159 CALL physics%new_statevector(this%cons,
conservative,mesh%NFACES)
162 CALL physics%new_statevector(this%pfluxes,
conservative,mesh%NFACES)
165 CALL this%Info(
" FLUXES---> fluxes type " // trim(this%GetName()))
168 CALL getattr(config,
"output/pstates", valwrite,0)
169 IF(valwrite.EQ.1)
THEN
171 key = trim(physics%pvarname(i)) //
"_pstates"
172 CALL setattr(io, trim(key), &
173 this%prim%data5d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,:,i))
176 CALL getattr(config,
"output/cstates", valwrite, 0)
177 IF(valwrite.EQ.1)
THEN
179 key = trim(physics%cvarname(i)) //
"_cstates"
180 CALL setattr(io, trim(key), &
181 this%cons%data5d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,:,i))
184 CALL getattr(config,
"output/pfluxes", valwrite, 0)
185 IF(valwrite.EQ.1)
THEN
187 key = trim(physics%pvarname(i)) //
"_pfluxes"
188 CALL setattr(io, trim(key), &
189 this%pfluxes%data5d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,:,i))
193 CALL getattr(config,
"output/wave_speeds", valwrite, 0)
194 IF(valwrite.EQ.1)
THEN
195 CALL setattr(io,
"minwav", &
196 this%minwav%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,1:mesh%NDIMS))
197 CALL setattr(io,
"maxwav", &
198 this%maxwav%data4d(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,1:mesh%NDIMS))
205 IF(
ASSOCIATED(iorec))
CALL setattr(io,
"reconstruction",iorec)
209 this%bxflux(:,:,:,:) = 0.
210 this%byflux(:,:,:,:) = 0.
211 this%bzflux(:,:,:,:) = 0.
212 this%bxfold(:,:,:,:) = 0.
213 this%byfold(:,:,:,:) = 0.
214 this%bzfold(:,:,:,:) = 0.
227 REAL,
DIMENSION(Physics%VNUM) :: bflux
228 INTEGER,
OPTIONAL :: comm
231 REAL,
DIMENSION(Physics%VNUM) :: bflux_local
232 INTEGER :: sender_rank(1),dest_ranks(1),rank0(1)
233 INTEGER :: dest_comm,union_comm
234 INTEGER :: world_group,dest_group,union_group,sender_group
238 INTENT(IN) :: direction
241 IF (
PRESENT(comm))
THEN
244 dest_comm = mpi_comm_null
247 SELECT CASE(direction)
250 IF (mesh%mycoords(1).EQ.0)
THEN
251 bflux_local(:) = sum(sum(
this%bxflux(mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,1,:),1),1)
256 bflux = sum(sum(
this%bxflux(mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,1,:),1),1)
260 IF (mesh%mycoords(1).EQ.mesh%dims(1)-1)
THEN
261 bflux_local(:) = sum(sum(
this%bxflux(mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,2,:),1),1)
266 bflux = sum(sum(
this%bxflux(mesh%JMIN:mesh%JMAX,mesh%KMIN:mesh%KMAX,2,:),1),1)
270 IF (mesh%mycoords(2).EQ.0)
THEN
271 bflux_local(:) = sum(sum(
this%byflux(mesh%KMIN:mesh%KMAX,mesh%IMIN:mesh%IMAX,1,:),1),1)
276 bflux = sum(sum(
this%byflux(mesh%KMIN:mesh%KMAX,mesh%IMIN:mesh%IMAX,1,:),1),1)
280 IF (mesh%mycoords(2).EQ.mesh%dims(2)-1)
THEN
281 bflux_local(:) = sum(sum(
this%byflux(mesh%KMIN:mesh%KMAX,mesh%IMIN:mesh%IMAX,2,:),1),1)
286 bflux = sum(sum(
this%byflux(mesh%KMIN:mesh%KMAX,mesh%IMIN:mesh%IMAX,2,:),1),1)
290 IF (mesh%mycoords(3).EQ.0)
THEN
291 bflux_local(:) = sum(sum(
this%bzflux(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,1,:),1),1)
296 bflux = sum(sum(
this%bzflux(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,1,:),1),1)
300 IF (mesh%mycoords(3).EQ.mesh%dims(3)-1)
THEN
301 bflux_local(:) = sum(sum(
this%bzflux(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,2,:),1),1)
306 bflux = sum(sum(
this%bzflux(mesh%IMIN:mesh%IMAX,mesh%JMIN:mesh%JMAX,2,:),1),1)
309 CALL this%Error(
"GetBoundaryFlux",
"wrong direction")
314 IF(dest_comm.EQ.mpi_comm_world)
THEN
316 mpi_sum,mpi_comm_world,ierror)
319 CALL mpi_comm_group(mpi_comm_world,world_group,ierror)
321 IF(dest_comm.NE.mpi_comm_null)
THEN
324 CALL mpi_comm_group(dest_comm,dest_group,ierror)
329 CALL mpi_group_incl(world_group,1,dest_ranks,dest_group,ierror)
333 IF (mesh%comm_boundaries(direction).NE.mpi_comm_null)
THEN
335 mpi_sum,0,mesh%comm_boundaries(direction),ierror)
340 rank0(1) = mesh%rank0_boundaries(direction)
341 CALL mpi_group_incl(world_group,1,rank0,sender_group,ierror)
343 CALL mpi_group_union(sender_group,dest_group,union_group,ierror)
345 CALL mpi_comm_create(mpi_comm_world,union_group,union_comm,ierror)
346 IF (union_comm.NE.mpi_comm_null)
THEN
349 CALL mpi_group_translate_ranks(sender_group,1,rank0,union_group,sender_rank,ierror)
350 IF (sender_rank(1).EQ.mpi_undefined) &
351 CALL this%Error(
"GetBoundaryFlux",
"sender rank undefined")
353 CALL mpi_bcast(bflux,physics%VNUM,
default_mpi_real,sender_rank(1),union_comm,ierror)
355 CALL mpi_comm_free(union_comm,ierror)
358 CALL mpi_group_free(union_group,ierror)
359 CALL mpi_group_free(sender_group,ierror)
360 CALL mpi_group_free(world_group,ierror)
361 CALL mpi_group_free(dest_group,ierror)
378 IF (this%Reconstruction%PrimRecon())
THEN
379 CALL this%Reconstruction%CalculateStates(mesh,physics,pvar,this%prim)
380 CALL physics%Convert2Conservative(this%prim,this%cons)
382 CALL this%Reconstruction%CalculateStates(mesh,physics,cvar,this%cons)
383 CALL physics%Convert2Primitive(this%cons,this%prim)
387 SELECT TYPE(phys => physics)
389 SELECT TYPE(prim => this%prim)
391 CALL phys%UpdateSoundSpeed(prim)
396 CALL physics%CalculateWaveSpeeds(mesh,this%prim%data5d,this%cons%data5d,this%minwav,this%maxwav)
405 IF (.NOT.this%Initialized()) &
406 CALL this%Error(
"CloseFluxes",
"not initialized")
407 DEALLOCATE(this%minwav,this%maxwav,this%prim,this%cons,this%pfluxes, &
408 this%bxflux,this%byflux,this%bzflux,this%bxfold,this%byfold,this%bzfold)
409 CALL this%Reconstruction%Finalize()
410 DEALLOCATE(this%Reconstruction)
Dictionary for generic data types.
type(logging_base), save this
base module for numerical flux functions
subroutine calculatefacedata(this, Mesh, Physics, pvar, cvar)
Calcualtes face data with reconstruction methods (e. g. limiters)
subroutine initfluxes(this, Mesh, Physics, config, IO, ftype, fname)
Initialize Fluxes.
subroutine finalize_base(this)
Destructor.
integer, parameter, public kt
real function, dimension(physics%vnum) getboundaryflux(this, Mesh, Physics, direction, comm)
Get fluxes at boundaries.
subroutine finalize(this)
Destructor of logging_base class.
integer, save default_mpi_real
default real type for MPI
base class for mesh arrays
derived class for compound of mesh arrays
integer, parameter east
named constant for eastern boundary
integer, parameter bottom
named constant for bottom boundary
integer, parameter south
named constant for southern boundary
integer, parameter top
named constant for top boundary
integer, parameter north
named constant for northern boundary
integer, parameter west
named constant for western boundary
physics module for 1D,2D and 3D non-isothermal Euler equations
constructor for physics class
constructor for reconstruction class
subroutine new_reconstruction(Reconstruction, Mesh, Physics, config, IO)