-
Notifications
You must be signed in to change notification settings - Fork 12
/
Copy pathex_2.f90
129 lines (84 loc) · 3.66 KB
/
ex_2.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
! gfortran -pedantic -cpp -O3 -flto ex_2.f90 -o ex_2
! ./ex_2
! Test for the monotone limiter: remap a stairstep profile
! between randomly perturbed grids. When MONO-LIMIT is
! selected, all methods should lead to monotone behaviour.
!
# include "../src/ppr_1d.f90"
program ex
use ppr_1d
implicit none
integer, parameter :: npos = 97 ! no. edge (old grid)
integer, parameter :: ntmp = 77 ! no. edge (new grid)
integer, parameter :: nvar = 1 ! no. variables to remap
integer, parameter :: ndof = 1 ! no. FV DoF per cell
integer :: ipos
!------------------------------ position of cell edges !
real(kind=dp) :: xpos(npos),xtmp(ntmp)
real(kind=dp) :: xmid
!-------------------------------- finite-volume arrays !
! Arrays represent a "block" of finite-volume tracers
! to remap. The 1st dim. is the no. of DoF per cell,
! NDOF=1 is a standard finite-volume scheme where the
! data is specified as cell means. NDOF>1 is reserved
! for future use with DG-style schemes. NVAR is the
! number of tracers to remap. Processing tracers in a
! batch is typically more efficient than one-by-one.
! The last dim. is the no. cells (layers) in the grid.
real(kind=dp) :: init(ndof,nvar,npos-1)
real(kind=dp) :: fdat(ndof,nvar,npos-1)
real(kind=dp) :: ftmp(ndof,nvar,ntmp-1)
!------------------------------ method data-structures !
type(rmap_work) :: work
type(rmap_opts) :: opts
type(rcon_ends) :: bc_l(nvar)
type(rcon_ends) :: bc_r(nvar)
!------------------------------ define a simple domain !
call linspace(0.d0,1.d0,npos,xpos)
call rndspace(0.d0,1.d0,ntmp,xtmp)
!------------------------------ setup some simple data !
do ipos = +1, npos-1
xmid = xpos(ipos+0) * 0.5d+0 &
& + xpos(ipos+1) * 0.5d+0
if (xmid .lt. 0.075d0) then
init(1,1,ipos) = +1.0d0
else &
if (xmid .lt. 0.80d0) then
init(1,1,ipos) = +2.0d0
else
init(1,1,ipos) = -0.5d0 * xmid ** 2
end if
end do
!------------------------------ specify method options !
opts%edge_meth = p5e_method ! 5th-order edge interp.
opts%cell_meth = pqm_method ! PQM method in cells
opts%cell_lims = mono_limit ! monotone limiter
!opts%wall_lims = weno_limit
!------------------------------ set BC.'s at endpoints !
bc_l%bcopt = bcon_loose ! "loose" = extrapolate
bc_r%bcopt = bcon_loose
!------------------------------ init. method workspace !
call work%init(npos,nvar,opts)
!------------------------------ re-map back-and-forth: !
fdat = init
do ipos = +1, +1000
!------------------------------ re-map from dat-to-tmp !
call rmap1d(npos,ntmp,nvar,ndof, &
& xpos,xtmp,fdat,ftmp, &
& bc_l,bc_r,work,opts)
!------------------------------ re-map from tmp-to-dat !
call rmap1d(ntmp,npos,nvar,ndof, &
& xtmp,xpos,ftmp,fdat, &
& bc_l,bc_r,work,opts)
end do
!------------------------------ clear method workspace !
call work%free()
!------------------------------ dump results to stdout !
print*,"Cell data: [INIT] [RMAP] "
do ipos = +1, npos-1
print *, init(1,1,ipos) &
& , fdat(1,1,ipos)
end do
print*,"Conservation defect := " &
& , sum(init) - sum(fdat)
end program