1 | |
---|
2 | ! |
---|
3 | ! This program may be freely redistributed under the |
---|
4 | ! condition that the copyright notices (including this |
---|
5 | ! entire header) are not removed, and no compensation |
---|
6 | ! is received through use of the software. Private, |
---|
7 | ! research, and institutional use is free. You may |
---|
8 | ! distribute modified versions of this code UNDER THE |
---|
9 | ! CONDITION THAT THIS CODE AND ANY MODIFICATIONS MADE |
---|
10 | ! TO IT IN THE SAME FILE REMAIN UNDER COPYRIGHT OF THE |
---|
11 | ! ORIGINAL AUTHOR, BOTH SOURCE AND OBJECT CODE ARE |
---|
12 | ! MADE FREELY AVAILABLE WITHOUT CHARGE, AND CLEAR |
---|
13 | ! NOTICE IS GIVEN OF THE MODIFICATIONS. Distribution |
---|
14 | ! of this code as part of a commercial system is |
---|
15 | ! permissible ONLY BY DIRECT ARRANGEMENT WITH THE |
---|
16 | ! AUTHOR. (If you are not directly supplying this |
---|
17 | ! code to a customer, and you are instead telling them |
---|
18 | ! how they can obtain it for free, then you are not |
---|
19 | ! required to make any arrangement with me.) |
---|
20 | ! |
---|
21 | ! Disclaimer: Neither I nor: Columbia University, the |
---|
22 | ! National Aeronautics and Space Administration, nor |
---|
23 | ! the Massachusetts Institute of Technology warrant |
---|
24 | ! or certify this code in any way whatsoever. This |
---|
25 | ! code is provided "as-is" to be used at your own risk. |
---|
26 | ! |
---|
27 | ! |
---|
28 | |
---|
29 | ! |
---|
30 | ! RCON1D.h90: conservative, polynomial reconstructions. |
---|
31 | ! |
---|
32 | ! Darren Engwirda |
---|
33 | ! 07-Sep-2016 |
---|
34 | ! de2363 [at] columbia [dot] edu |
---|
35 | ! |
---|
36 | ! |
---|
37 | |
---|
38 | subroutine rcon1d(npos,nvar,ndof,delx,fdat, & |
---|
39 | & bclo,bchi,fhat,work,opts, & |
---|
40 | & tCPU) |
---|
41 | |
---|
42 | ! |
---|
43 | ! NPOS no. edges over grid. |
---|
44 | ! NVAR no. state variables. |
---|
45 | ! NDOF no. degrees-of-freedom per grid-cell. |
---|
46 | ! DELX grid-cell spacing array. LENGTH(DELX) == +1 if |
---|
47 | ! spacing is uniform . |
---|
48 | ! FDAT grid-cell moments array. FDAT is an array with |
---|
49 | ! SIZE = NDOF-by-NVAR-by-NPOS-1 . |
---|
50 | ! BCLO boundary condition at lower endpoint. |
---|
51 | ! BCHI boundary condition at upper endpoint. |
---|
52 | ! FHAT grid-cell re-con. array. FHAT is an array with |
---|
53 | ! SIZE = MDOF-by-NVAR-by-NPOS-1 . |
---|
54 | ! WORK method work-space. See RCON-WORK for details . |
---|
55 | ! OPTS method parameters. See RCON-OPTS for details . |
---|
56 | ! TCPU method tcpu-timer. |
---|
57 | ! |
---|
58 | |
---|
59 | implicit none |
---|
60 | |
---|
61 | !------------------------------------------- arguments ! |
---|
62 | integer, intent( in) :: npos,nvar,ndof |
---|
63 | class(rcon_work), intent(inout):: work |
---|
64 | class(rcon_opts), intent(in) :: opts |
---|
65 | real*8 , intent( in) :: delx(:) |
---|
66 | real*8 , intent(out) :: fhat(:,:,:) |
---|
67 | real*8 , intent( in) :: fdat(:,:,:) |
---|
68 | type (rcon_ends), intent(in) :: bclo(:) |
---|
69 | type (rcon_ends), intent(in) :: bchi(:) |
---|
70 | type (rmap_tics), & |
---|
71 | & intent(inout) , optional :: tCPU |
---|
72 | |
---|
73 | !------------------------------------------- variables ! |
---|
74 | integer :: halo,ipos |
---|
75 | real*8 :: dmin,dmid |
---|
76 | |
---|
77 | # ifdef __PPR_TIMER__ |
---|
78 | integer(kind=8) :: ttic,ttoc,rate |
---|
79 | # endif |
---|
80 | |
---|
81 | if (ndof.lt.1) return |
---|
82 | if (npos.lt.2) return |
---|
83 | if (nvar.lt.1) return |
---|
84 | |
---|
85 | !-------------------------- compute min grid-tolerance ! |
---|
86 | |
---|
87 | dmid = delx(1) |
---|
88 | |
---|
89 | if (size(delx).gt.+1) then |
---|
90 | |
---|
91 | do ipos = 2, npos-1 |
---|
92 | dmid = & |
---|
93 | & dmid + delx (ipos) |
---|
94 | end do |
---|
95 | |
---|
96 | dmid = dmid /(npos-1) |
---|
97 | |
---|
98 | end if |
---|
99 | |
---|
100 | dmin = +1.0d-14 * dmid |
---|
101 | |
---|
102 | !-------------------------- compute edge values/slopes ! |
---|
103 | |
---|
104 | |
---|
105 | if ( (opts%cell_meth.eq.ppm_method) & |
---|
106 | & .or. (opts%cell_meth.eq.pqm_method) ) then |
---|
107 | |
---|
108 | select case (opts%edge_meth) |
---|
109 | case(p1e_method) |
---|
110 | !------------------------------------ 2nd-order method ! |
---|
111 | halo = +1 |
---|
112 | call p1e(npos,nvar,ndof, & |
---|
113 | & delx,fdat, & |
---|
114 | & bclo,bchi, & |
---|
115 | & work%edge_func, & |
---|
116 | & work%edge_dfdx, & |
---|
117 | & opts,dmin) |
---|
118 | |
---|
119 | case(p3e_method) |
---|
120 | !------------------------------------ 4th-order method ! |
---|
121 | halo = +2 |
---|
122 | call p3e(npos,nvar,ndof, & |
---|
123 | & delx,fdat, & |
---|
124 | & bclo,bchi, & |
---|
125 | & work%edge_func, & |
---|
126 | & work%edge_dfdx, & |
---|
127 | & opts,dmin) |
---|
128 | |
---|
129 | case(p5e_method) |
---|
130 | !------------------------------------ 6th-order method ! |
---|
131 | halo = +3 |
---|
132 | call p5e(npos,nvar,ndof, & |
---|
133 | & delx,fdat, & |
---|
134 | & bclo,bchi, & |
---|
135 | & work%edge_func, & |
---|
136 | & work%edge_dfdx, & |
---|
137 | & opts,dmin) |
---|
138 | |
---|
139 | end select |
---|
140 | |
---|
141 | end if |
---|
142 | |
---|
143 | |
---|
144 | !-------------------------- compute oscil. derivatives ! |
---|
145 | |
---|
146 | |
---|
147 | if (opts%cell_lims.eq.weno_limit) then |
---|
148 | |
---|
149 | call oscli(npos,nvar,ndof, & |
---|
150 | & delx,fdat, & |
---|
151 | & work%cell_oscl, & |
---|
152 | & dmin) |
---|
153 | |
---|
154 | end if |
---|
155 | |
---|
156 | |
---|
157 | !-------------------------- compute grid-cell profiles ! |
---|
158 | |
---|
159 | |
---|
160 | select case (opts%cell_meth) |
---|
161 | case(pcm_method) |
---|
162 | !------------------------------------ 1st-order method ! |
---|
163 | call pcm(npos,nvar,ndof, & |
---|
164 | & fdat,fhat) |
---|
165 | |
---|
166 | case(plm_method) |
---|
167 | !------------------------------------ 2nd-order method ! |
---|
168 | call plm(npos,nvar,ndof, & |
---|
169 | & delx,fdat,fhat, & |
---|
170 | & dmin,& |
---|
171 | & opts%cell_lims) |
---|
172 | |
---|
173 | case(ppm_method) |
---|
174 | !------------------------------------ 3rd-order method ! |
---|
175 | call ppm(npos,nvar,ndof, & |
---|
176 | & delx,fdat,fhat, & |
---|
177 | & work%edge_func, & |
---|
178 | & work%cell_oscl, & |
---|
179 | & dmin,& |
---|
180 | & opts%cell_lims, & |
---|
181 | & opts%wall_lims, & |
---|
182 | & halo ) |
---|
183 | |
---|
184 | case(pqm_method) |
---|
185 | !------------------------------------ 5th-order method ! |
---|
186 | call pqm(npos,nvar,ndof, & |
---|
187 | & delx,fdat,fhat, & |
---|
188 | & work%edge_func, & |
---|
189 | & work%edge_dfdx, & |
---|
190 | & work%cell_oscl, & |
---|
191 | & dmin,& |
---|
192 | & opts%cell_lims, & |
---|
193 | & opts%wall_lims, & |
---|
194 | & halo ) |
---|
195 | |
---|
196 | end select |
---|
197 | |
---|
198 | |
---|
199 | end subroutine |
---|
200 | |
---|
201 | |
---|
202 | |
---|