fosite.f90
Go to the documentation of this file.
1!#############################################################################
2!# #
3!# fosite - 3D hydrodynamical simulation program #
4!# program file: fosite.f90 #
5!# #
6!# Copyright (C) 2006-2024 #
7!# Tobias Illenseer <tillense@astrophysik.uni-kiel.de> #
8!# Manuel Jung <mjung@astrophysik.uni-kiel.de> #
9!# Björn Sperling <sperling@astrophysik.uni-kiel.de> #
10!# Jannes Klee <jklee@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!----------------------------------------------------------------------------!
40!----------------------------------------------------------------------------!
52 USE integration
53 USE common_dict
54#ifdef PARALLEL
55#ifdef HAVE_MPI_MOD
56 USE mpi
57#endif
58#endif
59 IMPLICIT NONE
60#ifdef PARALLEL
61#ifdef HAVE_MPIF_H
62 include 'mpif.h'
63#endif
64#endif
65 !--------------------------------------------------------------------------!
66 INTEGER, PRIVATE, PARAMETER :: maxlen = 500
67 INTEGER, PRIVATE, PARAMETER :: simtype = 1
68 CHARACTER(LEN=32), PRIVATE, PARAMETER :: simname = "fosite"
69
71 TYPE, EXTENDS(logging_base) :: fosite
74 TYPE(dict_typ),POINTER :: config => null()
75 TYPE(dict_typ),POINTER :: io => null()
76 CLASS(mesh_base),ALLOCATABLE :: mesh
77 CLASS(fluxes_base),ALLOCATABLE :: fluxes
78 CLASS(physics_base),ALLOCATABLE :: physics
79 CLASS(fileio_base),ALLOCATABLE :: datafile
80 CLASS(timedisc_base),ALLOCATABLE :: timedisc
81 CLASS(fileio_base),ALLOCATABLE :: logfile
82 CLASS(sources_list),ALLOCATABLE :: sources
85 INTEGER :: iter
86 LOGICAL :: aborted = .true.
87 DOUBLE PRECISION :: wall_time = 0.0
88 DOUBLE PRECISION :: log_time = 0.0
89 DOUBLE PRECISION :: start_time = 0.0
90 DOUBLE PRECISION :: end_time = 0.0
91 DOUBLE PRECISION :: run_time = 0.0
92 INTEGER :: start_count
93 CHARACTER(MAXLEN) :: buffer
96#ifdef PARALLEL
97 INTEGER :: ierror
98 REAL :: dt_all
100#endif
101
102 CONTAINS
105 PROCEDURE :: initfosite
106 PROCEDURE :: setup
107 PROCEDURE :: firststep
108 PROCEDURE :: step
109 PROCEDURE :: run
110 PROCEDURE :: printinfo
112 PROCEDURE :: printsummary
113 PROCEDURE :: computeruntime
114 PROCEDURE :: finalize
115 END TYPE fosite
116 !--------------------------------------------------------------------------!
117 PUBLIC :: fosite
120 !--------------------------------------------------------------------------!
121
122
123CONTAINS
124
125 SUBROUTINE initfosite(this)
126 IMPLICIT NONE
127 !--------------------------------------------------------------------------!
128 CLASS(fosite), INTENT(INOUT) :: this
129 !--------------------------------------------------------------------------!
130#ifdef PARALLEL
131 LOGICAL :: already_initialized = .false.
132#endif
133 !--------------------------------------------------------------------------!
134 IF (this%Initialized()) THEN
135 CALL this%Warning("InitFosite","already initialized, trying to reset")
136 ! finalize fosite without finalizing MPI
137 CALL this%Finalize(.false.)
138 END IF
139#ifdef PARALLEL
140 ! check whether MPI is already initialized
141 CALL mpi_initialized(already_initialized,this%ierror)
142 IF (.NOT.already_initialized) &
143 CALL mpi_init(this%ierror)
144#endif
145
146 CALL this%InitLogging(simtype,simname)
147
148 CALL initdict()
149 this%iter = -1
150 END SUBROUTINE initfosite
151
152 SUBROUTINE setup(this)
153 IMPLICIT NONE
154 !--------------------------------------------------------------------------!
155 CLASS(fosite), INTENT(INOUT) :: this
156 TYPE(dict_typ), POINTER :: dir, iodir, config_copy
157 TYPE(dict_typ), POINTER :: boundary
158 !--------------------------------------------------------------------------!
159 IF (.NOT.this%Initialized()) &
160 CALL this%Error("Setup","Sim is uninitialized")
161
162 CALL deletedict(this%IO)
163
164 IF (this%GetRank().EQ.0) THEN
165 ! print some information
166 WRITE(this%buffer, "(A)")&
167 "+---------------------------------------------------------+"
168 CALL this%Info(this%buffer)
169 WRITE(this%buffer, "(A1,A29,A28,A1)")&
170 "|",trim(simname),"","|"
171 CALL this%Info(this%buffer)
172 WRITE(this%buffer, "(A1,A35,A22,A1)")&
173 "|",trim(version),"","|"
174 CALL this%Info(this%buffer)
175 WRITE(this%buffer, "(A)")&
176 "| Solution of 3D advection problems |"
177 CALL this%Info(this%buffer)
178 WRITE(this%buffer, "(A)")&
179 "+---------------------------------------------------------+"
180 CALL this%Info(this%buffer)
181 CALL this%Info("Initializing simulation:")
182 END IF
183
184 CALL getattr(this%config, "mesh", dir)
185 IF(ASSOCIATED(dir)) THEN
186 NULLIFY(iodir)
187 CALL new_mesh(this%Mesh, dir, iodir)
188 IF(ASSOCIATED(iodir)) CALL setattr(this%IO, "mesh", iodir)
189 IF (this%Mesh%shear_dir.GT.0) THEN
190 ! add shearing box source term to the config
191 CALL setattr(this%config,"sources/shearing/stype",shearbox)
192 END IF
193 END IF
194
195 CALL getattr(this%config, "physics", dir)
196 IF(ASSOCIATED(dir)) THEN
197 NULLIFY(iodir)
198 CALL new_physics(this%Physics, this%Mesh, dir, iodir)
199 IF(ASSOCIATED(iodir)) CALL setattr(this%IO, "physics", iodir)
200 END IF
201
202 CALL getattr(this%config, "fluxes", dir)
203 IF(ASSOCIATED(dir)) THEN
204 NULLIFY(iodir)
205 CALL new_fluxes(this%Fluxes, this%Mesh, this%Physics, dir, iodir)
206 IF(ASSOCIATED(iodir)) CALL setattr(this%IO, "fluxes", iodir)
207 END IF
208
209 CALL getattr(this%config, "timedisc", dir)
210 IF(ASSOCIATED(dir)) THEN
211 NULLIFY(iodir)
212 CALL new_timedisc(this%Timedisc, this%Mesh,this%Physics,dir,iodir)
213 IF(ASSOCIATED(iodir)) CALL setattr(this%IO, "timedisc", iodir)
214 END IF
215
216 CALL getattr(this%config, "boundary", dir)
217 IF(.NOT.ASSOCIATED(dir)) THEN
218 boundary => dict("empty"/ 0)
219 CALL setattr(this%config, "boundary", boundary)
220 CALL getattr(this%config, "boundary", dir)
221 END IF
222 IF(ASSOCIATED(dir)) THEN
223 NULLIFY(iodir)
224 CALL new_boundary(this%Timedisc%Boundary, this%Mesh, this%Physics, dir,iodir)
225 IF(ASSOCIATED(iodir)) CALL setattr(this%IO, "boundary", iodir)
226 END IF
227
228 CALL getattr(this%config, "sources", dir)
229 IF(ASSOCIATED(dir)) THEN
230 IF(haschild(dir)) THEN
231 NULLIFY(iodir)
232 ALLOCATE(sources_list::this%Sources)
233 CALL this%Sources%InitSources(this%Mesh,this%Physics,this%Fluxes,dir,iodir)
234 IF(ASSOCIATED(iodir)) CALL setattr(this%IO, "sources", iodir)
235 END IF
236 END IF
237
238 CALL getattr(this%config, "datafile", dir)
239 IF(ASSOCIATED(dir)) THEN
240 CALL new_fileio(this%Datafile, this%Mesh, this%Physics, this%Timedisc,&
241 this%Sources, dir,this%IO)
242 END IF
243
244 CALL getattr(this%config, "logfile", dir)
245 IF(ASSOCIATED(dir)) THEN
246 CALL new_fileio(this%Logfile, this%Mesh, this%Physics, this%Timedisc,&
247 this%Sources,dir,this%IO)
248 END IF
249
250 CALL setattr(this%config, "version", trim(version))
251 CALL copydict(this%config, config_copy)
252
253 CALL setattr(this%IO, "config", config_copy)
254
255 END SUBROUTINE setup
256
257
258 SUBROUTINE firststep(this)
260 IMPLICIT NONE
261 !--------------------------------------------------------------------------!
262 CLASS(fosite), INTENT(INOUT) :: this
263 INTEGER :: datetime(8),dir
264 !--------------------------------------------------------------------------!
265 IF(.NOT.this%Initialized()) &
266 CALL this%Error("fosite::FirstStep","Sim is uninitialized")
267
268 ! this%iter = -1 initially (see InitFosite)
269 ! it is set to 0 here to make sure FirstStep is only called once
270 IF (this%iter.GT.-1) &
271 CALL this%Error("fosite::FirstStep","FirstStep should only be called once.")
272 this%iter = 0
273
274#ifdef PARALLEL
275 CALL mpi_barrier(mpi_comm_world,this%ierror)
276 this%start_time = mpi_wtime()
277 CALL mpi_allreduce(mpi_in_place,this%start_time,1,mpi_double_precision,mpi_max,&
278 mpi_comm_world,this%ierror)
279#else
280 CALL cpu_time(this%start_time)
281#endif
282
283 CALL system_clock(this%start_count)
284
285 ! make sure that the initial data is written to the log file
286 this%wall_time = this%start_time
287 this%log_time = this%wall_time
288
289 IF (this%GetRank().EQ.0) THEN
290 CALL this%Info( &
291 "===================================================================")
292 CALL this%Info("Starting calculation...")
293 CALL date_and_time(values = datetime)
294 WRITE(this%buffer,"(A,I2.2,A,I2.2,A,I4.4,A,I2.2,A,I2.2,A,I2.2)")&
295 "Time: ", &
296 datetime(3), ".", datetime(2), ".", datetime(1), " ",&
297 datetime(5), ":", datetime(6), ":", datetime(7)
298 CALL this%Info(this%buffer)
299 CALL this%Info( &
300 "step time n t min(dt) due to adj.")
301 CALL this%Info( &
302 "-------------------------------------------------------------------")
303 END IF
304
305 ! determine the background velocity if fargo advection type 1 is enabled
306 IF (this%Mesh%fargo%GetType().EQ.1) THEN
307 ! make sure there is valid data in the ghost cells
308 CALL this%Timedisc%Boundary%CenterBoundary(this%Mesh,this%Physics,&
309 0.0,this%Timedisc%pvar,this%Timedisc%cvar)
310 CALL this%Timedisc%CalcBackgroundVelocity(this%Mesh,this%Physics, &
311 this%Timedisc%pvar,this%Timedisc%cvar)
312 END IF
313
314 ! do a complete update of all data
315 CALL this%Timedisc%ComputeRHS(this%Mesh,this%Physics,this%Sources,this%Fluxes, &
316 this%Timedisc%time,0.0,this%Timedisc%pvar,this%Timedisc%cvar, &
317 check_all,this%Timedisc%rhs)
318
319 ! calculate timestep
320 this%Timedisc%dt = this%Timedisc%CalcTimestep(this%Mesh,this%Physics,this%Sources,&
321 this%Fluxes,this%Timedisc%time,this%Timedisc%dtcause)
322
323 SELECT TYPE(phys => this%Physics)
324 CLASS IS(physics_eulerisotherm)
325 IF(phys%csiso.GT.0.) THEN
326 IF(any(phys%bccsound%data1d(:).NE.phys%csiso)) THEN
327 CALL this%Error("FirstStep","isothermal sound speed set, but "&
328 // "arrays bccsound and/or fcsound have been overwritten.")
329 END IF
330 END IF
331 CLASS DEFAULT
332 ! do nothing
333 END SELECT
334
335
336 ! store old values
337 this%Timedisc%cold = this%Timedisc%cvar
338
339 ! store initial data
340 IF (this%Timedisc%time.EQ.0.0) THEN
341 CALL this%Datafile%WriteDataset(this%Mesh,this%Physics,this%Fluxes,&
342 this%Timedisc,this%config,this%IO)
343 CALL this%PrintInfo(0,0,0.0,0.0,0,this%Timedisc%n_adj)
344 END IF
345
346 IF(this%Timedisc%break) &
347 CALL this%Error("FirstStep","Initial data invalid!")
348 END SUBROUTINE firststep
349
350
351 SUBROUTINE run(this)
352 IMPLICIT NONE
353 !--------------------------------------------------------------------------!
354 CLASS(fosite), INTENT(INOUT) :: this
355 !--------------------------------------------------------------------------!
356 CHARACTER(LEN=32) :: str
357 !--------------------------------------------------------------------------!
358 ! main loop
359 DO WHILE((this%Timedisc%maxiter.LE.0).OR.(this%iter.LE.this%Timedisc%maxiter))
360 IF(this%Step()) EXIT
361 END DO
362 IF (this%iter.GE.this%Timedisc%maxiter) &
363 CALL this%Warning("fosite::Run","number of iterations exceeds limit")
364 END SUBROUTINE run
365
366
367 FUNCTION step(this) RESULT(break)
368 IMPLICIT NONE
369 !--------------------------------------------------------------------------!
370 CLASS(fosite), INTENT(INOUT) :: this
371 LOGICAL :: break
372#ifdef PARALLEL
373 REAL,DIMENSION(2) :: dt_buf
374#endif
375! TYPE(Gravity_TYP), POINTER :: pmass
376 !--------------------------------------------------------------------------!
377 break = .false.
378
379 IF (this%iter.EQ.-1) &
380 CALL this%FirstStep()
381
382 ! finish simulation if stop time is reached
383 IF (abs(this%Timedisc%stoptime-this%Timedisc%time)&
384 .LE.1.0e-05*this%Timedisc%stoptime) THEN
385 break = .true.
386 this%aborted = .false.
387 RETURN
388 END IF
389
390#ifdef PARALLEL
391 ! In Fortran MPI_MINLOC is only able to have two values of the same kind!
392 dt_buf(1) = this%Timedisc%dt
393 dt_buf(2) = this%Timedisc%dtcause
394 CALL mpi_allreduce(mpi_in_place,dt_buf,1,default_mpi_2real,mpi_minloc,&
395 this%Mesh%comm_cart,this%ierror)
396 this%Timedisc%dt = dt_buf(1)
397 this%Timedisc%dtcause = dt_buf(2)
398
399 this%run_time = mpi_wtime() - this%log_time
400 CALL mpi_allreduce(this%run_time,this%wall_time,1,mpi_double_precision,&
401 mpi_min,this%Mesh%comm_cart,this%ierror)
402#else
403 CALL cpu_time(this%wall_time)
404 this%wall_time = this%wall_time - this%log_time
405#endif
406
407 ! adjust timestep for output and calculate the wall clock time
408 CALL this%Datafile%AdjustTimestep(this%Timedisc%time,&
409 this%Timedisc%dt,this%Timedisc%dtcause)
410
411#ifdef PARALLEL
412 ! In Fortran MPI_MINLOC is only able to have two values of the same kind!
413 dt_buf(1) = this%Timedisc%dt
414 dt_buf(2) = this%Timedisc%dtcause
415 CALL mpi_allreduce(mpi_in_place,dt_buf,1,default_mpi_2real,mpi_minloc,&
416 this%Mesh%comm_cart,this%ierror)
417 this%Timedisc%dt = dt_buf(1)
418 this%Timedisc%dtcause = dt_buf(2)
419#endif
420
421 ! advance the solution in time
422 CALL this%Timedisc%IntegrationStep(this%Mesh,this%Physics,this%Sources, &
423 this%Fluxes,this%iter,this%config,this%IO)
424
425 ! write output to data file
426 IF ((abs(this%Datafile%time-this%Timedisc%time)&
427 .LE.1.0e-5*this%Datafile%time).OR.&
428 this%Timedisc%break.OR.(this%Datafile%walltime.LE.this%wall_time)) THEN
429 CALL this%Datafile%WriteDataset(this%Mesh,this%Physics,this%Fluxes,&
430 this%Timedisc,this%config,this%IO)
431 CALL this%PrintInfo(this%Datafile%step-1, this%iter,this%Timedisc%time,&
432 this%Timedisc%dtmin,this%Timedisc%dtmincause,this%Timedisc%n_adj)
433
434 ! Stop program if a break is requested from SolveODE
435 IF(this%Timedisc%break.AND.this%GetRank().EQ.0) THEN
436 CALL this%Warning("SolveODE", "Time step too small, aborting.",0)
437 break = .true.
438 RETURN
439 END IF
440
441 ! only give one additional output at before walltime
442 IF (this%Datafile%walltime.LE.this%wall_time) THEN
443 this%Datafile%walltime = huge(1.0)
444 END IF
445
446 ! reset dt_min,dtmincause and n_adj
447 this%Timedisc%dtmin = this%Timedisc%stoptime
448 this%Timedisc%dtmincause = -99
449 this%Timedisc%n_adj = 0
450 !IF(GetRank(this).EQ.0) &
451 ! WRITE(*,"(A,ES10.4,A,ES10.4)") "dtmean: ", this%Timedisc%dtmean, " +- ",&
452 ! SQRT(this%Timedisc%dtstddev/(this%Timedisc%dtaccept-1))/this%Timedisc%dtmean
453 this%Timedisc%dtmean = 0.
454 this%Timedisc%dtstddev = 0.
455 this%Timedisc%dtaccept = 0
456 ! reset max error of cvar
457 IF(this%Timedisc%write_error) &
458 this%Timedisc%cerr_max%data1d(:) = 0.
459 END IF
460
461 ! calculate next timestep
462 this%Timedisc%dt = this%Timedisc%CalcTimestep(this%Mesh,this%Physics,this%Sources, &
463 this%Fluxes,this%Timedisc%time,this%Timedisc%dtcause)
464 END FUNCTION step
465
466 SUBROUTINE finalize(this,mpifinalize_)
467 IMPLICIT NONE
468 !--------------------------------------------------------------------------!
469 CLASS(fosite), INTENT(INOUT) :: this
470 LOGICAL,OPTIONAL, INTENT(IN) :: mpifinalize_
471 !--------------------------------------------------------------------------!
472#ifdef PARALLEL
473 LOGICAL :: mpifinalize
474#endif
475 !--------------------------------------------------------------------------!
476 CALL this%ComputeRunTime()
477 CALL this%PrintBoundaryFluxes()
478 CALL this%PrintSummary()
479
480! CALL this%Datafile%Finalize()
481 DEALLOCATE(this%Datafile) ! should call finalizer automatically
482 IF (ALLOCATED(this%Logfile)) THEN
483! CALL this%Logfile%Finalize()
484 DEALLOCATE(this%Logfile) ! should call finalizer automatically
485 END IF
486 CALL this%Timedisc%Finalize()
487 DEALLOCATE(this%Timedisc)
488
489 IF (ALLOCATED(this%Sources)) DEALLOCATE(this%Sources)
490
491 CALL this%Fluxes%Finalize()
492 DEALLOCATE(this%Fluxes)
493
494 CALL this%Physics%Finalize()
495 DEALLOCATE(this%Physics)
496
497 CALL this%Mesh%Finalize()
498 DEALLOCATE(this%Mesh)
499
500 CALL deletedict(this%IO)
501 CALL deletedict(this%config)
502 CALL closedict()
503
504#ifdef PARALLEL
505 IF(PRESENT(mpifinalize_)) THEN
506 mpifinalize = mpifinalize_
507 ELSE
508 mpifinalize = .true.
509 END IF
510 IF(mpifinalize) &
511 CALL mpi_finalize(this%ierror)
512#endif
513 END SUBROUTINE finalize
514
515 SUBROUTINE computeruntime(this)
516 IMPLICIT NONE
517 !--------------------------------------------------------------------------!
518 CLASS(fosite), INTENT(INOUT) :: this
519 !--------------------------------------------------------------------------!
520
521#ifdef PARALLEL
522 CALL mpi_barrier(mpi_comm_world,this%ierror)
523 this%end_time = mpi_wtime()
524 CALL mpi_allreduce(mpi_in_place,this%end_time,1,mpi_double_precision,mpi_min,&
525 mpi_comm_world,this%ierror)
526#else
527 CALL cpu_time(this%end_time)
528#endif
529 this%run_time = this%end_time - this%start_time
530 END SUBROUTINE computeruntime
531
532
533 SUBROUTINE printinfo(this,step,i,t,d,dc,na)
534 IMPLICIT NONE
535 !------------------------------------------------------------------------!
536 CLASS(fosite), INTENT(INOUT) :: this
537 INTEGER :: step, i, dc, na
538 CHARACTER(LEN=9) :: dtcause
539 REAL :: t, d
540 !------------------------------------------------------------------------!
541 INTEGER :: drt, c, c_rate, c_max
542 !------------------------------------------------------------------------!
543 INTENT(IN) :: i,t,d
544 !------------------------------------------------------------------------!
545 IF (this%GetRank().EQ.0) THEN
546
547 SELECT CASE (dc)
548 ! positive values represent source terms
549 CASE(1:)
550 WRITE(dtcause, "(A,I2.2,A)") " S",dc," "
551 CASE(dtcause_cfl)
552 WRITE(dtcause, "(A)") " cfl "
553 CASE(dtcause_erradj)
554 WRITE(dtcause, "(A)") " err_adj "
555 CASE(dtcause_smallerr)
556 WRITE(dtcause, "(A)") " err "
557 CASE(dtcause_fileio)
558 ! output by this reason is suppressed by default
559 WRITE(dtcause, "(A)") " fileio "
560 CASE DEFAULT
561 WRITE(dtcause, "(A,I3.2,A)") " ?", dc, "? "
562 END SELECT
563
564 CALL system_clock(c,c_rate, c_max)
565 ! overflow every 24days => assumption: max one per simulation
566 drt = (c - this%start_count)/c_rate
567 IF (drt .LT. 0) THEN
568 drt = c_max/c_rate + drt
569 END IF
570 WRITE(this%buffer,"(I4.4,A,I2.2,A,I2.2,A,I2.2,A,I8,A,ES11.3,A,ES11.3,A,I5)")&
571 step, " ", drt/3600, ":", mod(drt,3600)/60, ":", mod(drt,60),&
572 " ", i, " ", t, " ", d, dtcause, na
573 CALL this%Info(this%buffer)
574 END IF
575 END SUBROUTINE printinfo
576
577 SUBROUTINE printboundaryfluxes(this)
578 IMPLICIT NONE
579 !------------------------------------------------------------------------!
580 CLASS(fosite), INTENT(INOUT) :: this
581 INTEGER :: k
582 REAL, DIMENSION(this%Physics%VNUM,6) :: bflux
583 !--------------------------------------------------------------------------!
584 DO k=1,6
585 bflux(:,k) = this%Fluxes%GetBoundaryFlux(this%Mesh,this%Physics,k)
586 END DO
587 IF (this%GetRank().EQ.0) THEN
588 CALL this%Info("-------------------------------------------------------------------")
589 CALL this%Info("total boundary fluxes:")
590 SELECT CASE(this%Physics%VDIM)
591 CASE(1)
592 IF(this%Mesh%INUM.GT.1) THEN
593 CALL this%Info(" west east")
594 DO k=1,this%Physics%VNUM
595 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
596 bflux(k,west), bflux(k,east)
597 CALL this%Info(this%buffer)
598 END DO
599 ELSEIF(this%Mesh%JNUM.GT.1) THEN
600 CALL this%Info(" south north")
601 DO k=1,this%Physics%VNUM
602 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
603 bflux(k,south), bflux(k,north)
604 CALL this%Info(this%buffer)
605 END DO
606 ELSEIF(this%Mesh%KNUM.GT.1) THEN
607 CALL this%Info(" bottom top")
608 DO k=1,this%Physics%VNUM
609 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
610 bflux(k,bottom), bflux(k,top)
611 CALL this%Info(this%buffer)
612 END DO
613 END IF
614 CASE(2)
615 IF(this%Mesh%INUM.EQ.1) THEN
616 CALL this%Info(" south north")
617 DO k=1,this%Physics%VNUM
618 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
619 bflux(k,south), bflux(k,north)
620 CALL this%Info(this%buffer)
621 END DO
622 CALL this%Info(" bottom top")
623 DO k=1,this%Physics%VNUM
624 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
625 bflux(k,bottom), bflux(k,top)
626 CALL this%Info(this%buffer)
627 END DO
628 ELSEIF(this%Mesh%JNUM.EQ.1) THEN
629 CALL this%Info(" west east")
630 DO k=1,this%Physics%VNUM
631 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
632 bflux(k,west), bflux(k,east)
633 CALL this%Info(this%buffer)
634 END DO
635 CALL this%Info(" bottom top")
636 DO k=1,this%Physics%VNUM
637 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
638 bflux(k,bottom), bflux(k,top)
639 CALL this%Info(this%buffer)
640 END DO
641 ELSEIF(this%Mesh%KNUM.EQ.1) THEN
642 CALL this%Info(" west east")
643 DO k=1,this%Physics%VNUM
644 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
645 bflux(k,west), bflux(k,east)
646 CALL this%Info(this%buffer)
647 END DO
648 CALL this%Info(" south north")
649 DO k=1,this%Physics%VNUM
650 WRITE(this%buffer,"(T2,A,T21,2(ES12.3))")trim(this%Physics%cvarname(k)), &
651 bflux(k,south), bflux(k,north)
652 CALL this%Info(this%buffer)
653 END DO
654 END IF
655 CASE(3)
656 CALL this%Info(" west east south")
657 DO k=1,this%Physics%VNUM
658 WRITE(this%buffer,"(T2,A,T21,3(ES12.3))")trim(this%Physics%cvarname(k)), &
659 bflux(k,west), bflux(k,east), bflux(k,south)
660 CALL this%Info(this%buffer)
661 END DO
662 CALL this%Info(" north bottom top")
663 DO k=1,this%Physics%VNUM
664 WRITE(this%buffer,"(T2,A,T21,3(ES12.3))")trim(this%Physics%cvarname(k)), &
665 bflux(k,north), bflux(k,bottom), bflux(k,top)
666 CALL this%Info(this%buffer)
667 END DO
668 END SELECT
669 END IF
670 END SUBROUTINE printboundaryfluxes
671
672 SUBROUTINE printsummary(this)
673 IMPLICIT NONE
674 !------------------------------------------------------------------------!
675 CLASS(fosite), INTENT(INOUT) :: this
676 !--------------------------------------------------------------------------!
677 IF (this%GetRank().EQ.0) THEN
678 CALL this%Info("===================================================================")
679 IF (this%aborted) THEN
680 CALL this%Info("time integration aborted due to errors!")
681 ELSE IF ((this%Timedisc%maxiter.LE.0).OR.(this%iter.LT.this%Timedisc%maxiter)) THEN
682 CALL this%Info("calculation finished correctly.")
683 ELSE
684 CALL this%Info("too many iterations, aborting!")
685 END IF
686 WRITE(this%buffer,"(A,F10.2,A)")" main loop runtime: ", this%run_time, " sec."
687 CALL this%Info(this%buffer)
688 END IF
689 END SUBROUTINE printsummary
690
691END MODULE fosite_mod
Generic boundary module.
subroutine, public new_boundary(Boundary, Mesh, Physics, config, IO)
Dictionary for generic data types.
Definition: common_dict.f90:61
recursive subroutine, public deletedict(root)
Delete the dictionary 'root' and all subnodes.
logical function, public haschild(root)
Check if the node 'root' has one or more children.
recursive subroutine, public copydict(root, outdir)
Copy complete Dictionary.
type(dict_typ) function, pointer, public dict(n1, n2, n3, n4, n5, n6, n7, n8, n9, n10, n11, n12, n13, n14, n15, n16, n17, n18, n19, n20)
Construct a new dictionary from several key/value pairs. Together with the Assign subroutine and over...
subroutine, public initdict()
subroutine, public closedict()
type(logging_base), save this
constructor for fileio class
subroutine new_fileio(Fileio, Mesh, Physics, Timedisc, Sources, config, IO)
constructor for FileIO class
constructor for fluxes class
subroutine, public new_fluxes(Fluxes, Mesh, Physics, config, IO)
subroutine, private run(this)
Definition: fosite.f90:352
logical function, private step(this)
Definition: fosite.f90:368
subroutine, private printinfo(this, step, i, t, d, dc, na)
Definition: fosite.f90:534
subroutine, private computeruntime(this)
Definition: fosite.f90:516
subroutine, private printboundaryfluxes(this)
Definition: fosite.f90:578
integer, parameter, private simtype
Definition: fosite.f90:67
character(len=32), parameter, private simname
Definition: fosite.f90:68
subroutine, private finalize(this, mpifinalize_)
Definition: fosite.f90:467
subroutine, private initfosite(this)
Definition: fosite.f90:126
subroutine, private printsummary(this)
Definition: fosite.f90:673
subroutine, private firststep(this)
Definition: fosite.f90:259
integer, parameter, private maxlen
Definition: fosite.f90:66
subroutine, private setup(this)
Definition: fosite.f90:153
defines properties of a 3D spherical mesh
Numerical integration.
Definition: integration.f90:33
Basic fosite module.
integer, save default_mpi_2real
default 2real type for MPI
constructor for mesh class
subroutine new_mesh(Mesh, config, IO)
constructor for physics class
subroutine new_physics(Physics, Mesh, config, IO)
allocate and initialize new physics class
constructor for reconstruction class
module to manage list of source terms
integer, parameter, public shearbox
constructor for timedisc class
subroutine new_timedisc(Timedisc, Mesh, Physics, config, IO)
main fosite class
Definition: fosite.f90:71
common data structure
container class to manage the list of source terms