RTI.f90
Go to the documentation of this file.
1 !#############################################################################
2 !# #
3 !# fosite - 3D hydrodynamical simulation program #
4 !# module: RTI.f90 #
5 !# #
6 !# Copyright (C) 2008-2018 #
7 !# Björn Sperling <sperling@astrophysik.uni-kiel.de> #
8 !# Tobias Illenseer <tillense@astrophysik.uni-kiel.de> #
9 !# Jannes Klee <jklee@astrophysik.uni-kiel.de> #
10 !# #
11 !# This program is free software; you can redistribute it and/or modify #
12 !# it under the terms of the GNU General Public License as published by #
13 !# the Free Software Foundation; either version 2 of the License, or (at #
14 !# your option) any later version. #
15 !# #
16 !# This program is distributed in the hope that it will be useful, but #
17 !# WITHOUT ANY WARRANTY; without even the implied warranty of #
18 !# MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE, GOOD TITLE or #
19 !# NON INFRINGEMENT. See the GNU General Public License for more #
20 !# details. #
21 !# #
22 !# You should have received a copy of the GNU General Public License #
23 !# along with this program; if not, write to the Free Software #
24 !# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. #
25 !# #
26 !#############################################################################
27 !----------------------------------------------------------------------------!
48 !----------------------------------------------------------------------------!
49 PROGRAM rti
50  USE fosite_mod
51 #include "tap.h"
52  IMPLICIT NONE
53  !--------------------------------------------------------------------------!
54  ! simulation parameters
55  REAL, PARAMETER :: tsim = 10.0 ! simulation time !
56  REAL, PARAMETER :: dynvis = 0.0 ! dynamic viscosity constant !
57  REAL, PARAMETER :: bulkvis = 0.0 ! bulk viscosity constant !
58 ! REAL, PARAMETER :: DYNVIS = 1.0E-4
59 ! REAL, PARAMETER :: BULKVIS = -6.67E-5
60  ! initial condition (SI units)
61  REAL, PARAMETER :: rho0 = 2.0 ! density: upper region !
62  REAL, PARAMETER :: rho1 = 1.0 ! density: lower region !
63  REAL, PARAMETER :: yacc = 0.2 ! grav. acceleration !
64  REAL, PARAMETER :: p0 = 1.2 ! pressure at the top !
65  REAL, PARAMETER :: a0 = 0.02 ! amplitude of init. disturbance !
66  ! mesh settings
67  INTEGER, PARAMETER :: xres = 50 ! resolution in x !
68  INTEGER, PARAMETER :: yres = 100 ! resolution in y !
69  INTEGER, PARAMETER :: zres = 1 ! z-resolution
70  REAL, PARAMETER :: width = 1.0 ! width of comp. domain !
71  REAL, PARAMETER :: height = 2.0 ! height of comp. domaina !
72  ! output parameters
73  INTEGER, PARAMETER :: onum = 10 ! number of output data sets !
74  CHARACTER(LEN=256), PARAMETER &
75  :: odir = './' ! output data dir !
76  CHARACTER(LEN=256), PARAMETER &
77  :: ofname = 'RTI' ! output data file name !
78  !--------------------------------------------------------------------------!
79  CLASS(fosite), ALLOCATABLE :: sim
80  !--------------------------------------------------------------------------!
81 
82 ! TAP_PLAN(1)
83 
84  ALLOCATE(sim)
85 
86  CALL sim%InitFosite()
87  CALL makeconfig(sim, sim%config)
88  CALL sim%Setup()
89  CALL initdata(sim%Mesh, sim%Physics, sim%Timedisc%pvar, sim%Timedisc%cvar)
90  CALL sim%Run()
91  CALL sim%Finalize()
92  DEALLOCATE(sim)
93 
94 ! TAP_CHECK(.TRUE.,"Simulation finished")
95 ! TAP_DONE
96 
97  CONTAINS
98 
100  SUBROUTINE makeconfig(Sim, config)
101  IMPLICIT NONE
102  !------------------------------------------------------------------------!
103  CLASS(fosite) :: Sim
104  TYPE(Dict_TYP), POINTER :: config
105  !------------------------------------------------------------------------!
106  TYPE(Dict_TYP), POINTER :: mesh, physics, boundary, datafile, &
107  sources, timedisc, fluxes, caccel, vis
108  !------------------------------------------------------------------------!
109  INTENT(INOUT) :: sim
110  !------------------------------------------------------------------------!
111  ! mesh settings
112  mesh => dict( &
113  "meshtype" / midpoint, & ! use midpoint rule !
114  "geometry" / cartesian, & ! cartesian grid !
115  "inum" / xres, & ! resolution in x-direction !
116  "jnum" / yres, & ! resolution in y-direction !
117  "knum" / zres, & ! resolution in z-direction !
118  "xmin" / 0., & ! minimum value in x-dir. !
119  "xmax" / width, & ! maximum value in x-dir. !
120  "ymin" / 0., & ! minimum value in y-dir. !
121  "ymax" / height, & ! maximum value in y-dir. !
122  "zmin" / 0., & ! minimum value in z-dir. !
123  "zmax" / 0.) ! maximum value in z-dir. !
124 
125  ! physics settings
126  physics => dict( &
127  "problem" / euler, & ! standard 2D hydrodynamics !
128  "gamma" / 1.4) ! ratio of specific heats !
129 
130  ! flux calculation and reconstruction method
131  fluxes => dict( &
132  "fluxtype" / kt, & ! use Kurganov-Tadmor flux !
133  "order" / linear, & ! linear reconstruction !
134  "variables" / conservative, & ! vars. for reconstruction !
135  "limiter" / monocent, & ! type of the limiter !
136  "theta" / 1.2) ! optional param. for MC !
137 
138  ! boundary conditions
139  boundary => dict( &
140  "western" / reflecting, & ! reflecting boundary cond. !
141  "eastern" / reflecting, & ! reflecting boundary cond. !
142  "southern" / reflecting, & ! reflecting boundary cond. !
143  "northern" / reflecting, & ! reflecting boundary cond. !
144  "bottomer" / reflecting, & ! reflecting boundary cond. !
145  "topper" / reflecting) ! reflecting boundary cond. !
146 
147  ! c-acceleration term
148  caccel => dict( &
149  "stype" / c_accel, & ! source 1: acceleration !
150  "yaccel" / (-yacc)) ! constant acc.in y-dir. !
151 
152 
153  ! viscosity source term
154  vis => dict( &
155  "stype" / viscosity, & ! viscosity source !
156  "vismodel" / molecular, & ! visc. model: molecular !
157  "dynconst" / dynvis, & ! const. dynamic viscosity !
158  "bulkconst" / bulkvis) ! const. bulk viscosity !
159 
160  ! collect sources in dictionary
161  sources => dict( &
162 ! "vis" / vis, & ! incl. visc. from above !
163  "caccel" / caccel) ! incl. accel. from above !
164 
165  IF ((dynvis.GT.tiny(dynvis)).OR.(bulkvis.GT.tiny(bulkvis))) &
166  CALL setattr(sources, "vis", vis)
167 
168  ! time discretization settings
169  timedisc => dict( &
170  "method" / modified_euler, &! use modified euler !
171  "order" / 3, & ! third order accuracy !
172  "cfl" / 0.4, & ! courant-number !
173  "stoptime" / tsim, & ! simulation stop-time !
174  "dtlimit" / 1.0e-4, & ! smallest allowed timestep !
175  "maxiter" / 100000) ! max. iters before abort !
176 
177  ! initialize data input/output
178  datafile => dict( &
179  "fileformat" / vtk, &
180  "filename" / (trim(odir) // trim(ofname)), &
181  "count" / onum)
182 
183  ! collect all above dicts in the configuration dict
184  config => dict( &
185  "mesh" / mesh, & ! all mesh-settings !
186  "physics" / physics, & ! all physics settings !
187  "boundary" / boundary, & ! all bounary settings !
188  "fluxes" / fluxes, & ! all fluxes settings !
189  "sources" / sources, & ! all sources !
190  "timedisc" / timedisc, & ! all timedisc settings !
191  "datafile" / datafile) ! all input/output settings !
192  END SUBROUTINE makeconfig
193 
194 
196  SUBROUTINE initdata(Mesh,Physics,pvar,cvar)
197  IMPLICIT NONE
198  !------------------------------------------------------------------------!
199  CLASS(mesh_base), INTENT(IN) :: Mesh
200  CLASS(physics_base), INTENT(IN) :: Physics
201  CLASS(marray_compound), POINTER, INTENT(INOUT) :: pvar,cvar
202  !------------------------------------------------------------------------!
203  ! Local variable declaration
204  REAL, DIMENSION(Mesh%IGMIN:Mesh%IGMAX,Mesh%JGMIN:Mesh%JGMAX,Mesh%KGMIN:Mesh%KGMAX) &
205  :: y0
206  !------------------------------------------------------------------------!
207  ! this marks the line between the two fluids
208  y0(:,:,:) = 0.5*mesh%ymax + a0*cos(2*pi*mesh%bcenter(:,:,:,1)/mesh%xmax)
209 
210  ! initial hydrostatic stratification
211  SELECT TYPE(p => pvar)
212  TYPE IS(statevector_euler)
213  WHERE (mesh%bcenter(:,:,:,2).GT.y0(:,:,:))
214  ! upper fluid
215  p%density%data3d(:,:,:) = rho0
216  p%pressure%data3d(:,:,:) = p0 + yacc * rho0 * (mesh%ymax-mesh%bcenter(:,:,:,2))
217  ELSEWHERE
218  ! lower fluid
219  p%density%data3d(:,:,:) = rho1
220  p%pressure%data3d(:,:,:) = p0 + yacc * (rho1 * (y0(:,:,:)-mesh%bcenter(:,:,:,2)) &
221  + rho0 * (mesh%ymax-y0(:,:,:)))
222  END WHERE
223  ! velocity vanishes everywhere
224  p%velocity%data1d(:) = 0.
225  CLASS DEFAULT
226  CALL physics%Error("shear::InitData","only non-isothermal HD supported")
227  END SELECT
228 
229  CALL physics%Convert2Conservative(pvar,cvar)
230  CALL mesh%Info(" DATA-----> initial condition: " // &
231  "Rayleigh–Taylor instability")
232  END SUBROUTINE initdata
233 END PROGRAM rti
subroutine initdata(Mesh, Physics, Fluxes, Timedisc)
Definition: bondi2d.f90:274
main fosite class
Definition: fosite.f90:71
subroutine makeconfig(Sim, config)
Definition: bondi2d.f90:165
program rti
Definition: RTI.f90:49