From 99cc1ca2abcf3e52ecfcf6bc66ee051cf034bbce Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Sun, 10 May 2015 08:09:43 -0300 Subject: [PATCH 1/8] INTERPOLATIONS: Make sure MC, VL, and SB limiters do not overshoot. Signed-off-by: Grzegorz Kowal --- src/interpolations.F90 | 18 ++++++++++++++---- 1 file changed, 14 insertions(+), 4 deletions(-) diff --git a/src/interpolations.F90 b/src/interpolations.F90 index daf671d..5970251 100644 --- a/src/interpolations.F90 +++ b/src/interpolations.F90 @@ -2483,10 +2483,15 @@ module interpolations ! real(kind=8), intent(in) :: x, a, b real(kind=8) :: c + +! local variables +! + real(kind=8) :: y ! !------------------------------------------------------------------------------- ! - c = (sign(x, a) + sign(x, b)) * min(abs(a), abs(b), 2.5d-01 * abs(a + b)) + y = x - eps + c = (sign(y, a) + sign(y, b)) * min(abs(a), abs(b), 2.5d-01 * abs(a + b)) !------------------------------------------------------------------------------- ! @@ -2517,11 +2522,16 @@ module interpolations ! real(kind=8), intent(in) :: x, a, b real(kind=8) :: c + +! local variables +! + real(kind=8) :: y ! !------------------------------------------------------------------------------- ! - c = 0.5d+00 * (sign(x, a) + sign(x, b)) & - * max(min(2.0d+00 * abs(a), abs(b)), min(abs(a), 2.0d+00 * abs(b))) + y = x - eps + c = (sign(y, a) + sign(y, b)) & + * max(min(abs(a), 0.5d+00 * abs(b)), min(0.5d+00 * abs(a), abs(b))) !------------------------------------------------------------------------------- ! @@ -2557,7 +2567,7 @@ module interpolations ! c = a * b if (c > 0.0d+00) then - c = 2.0d+00 * x * c / (a + b) + c = 2.0d+00 * (x - eps) * c / (a + b) else c = 0.0d+00 end if From 9c7962915d5e667b8c73a1c6eff12356e42ca56c Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Sun, 10 May 2015 08:33:26 -0300 Subject: [PATCH 2/8] INTERPOLATIONS: Apply TVD reconstruction to clipped extrema. Signed-off-by: Grzegorz Kowal --- src/interpolations.F90 | 33 ++++++++++++++++++++++++++++----- 1 file changed, 28 insertions(+), 5 deletions(-) diff --git a/src/interpolations.F90 b/src/interpolations.F90 index 5970251..90ca9ee 100644 --- a/src/interpolations.F90 +++ b/src/interpolations.F90 @@ -2729,8 +2729,9 @@ module interpolations ! local variables ! - integer :: i, im1, ip1 + integer :: i, im1, ip1, ip2 real(kind=8) :: fmn, fmx + real(kind=8) :: dfl, dfr, df ! !------------------------------------------------------------------------------ ! @@ -2758,10 +2759,19 @@ module interpolations ! if (fl(i) < fmn .or. fl(i) > fmx) then +! calculate the left and right derivatives +! + dfl = f(i ) - f(im1) + dfr = f(ip1) - f(i ) + +! get the limited slope +! + df = limiter_minmod(0.5d+00, dfl, dfr) + ! calculate new states ! - fl(i ) = f(i ) - fr(im1) = f(i ) + fl(i ) = f(i ) + df + fr(im1) = f(i ) - df end if @@ -2769,10 +2779,23 @@ module interpolations ! if (fr(i) < fmn .or. fr(i) > fmx) then +! calculate the missing index +! + ip2 = min(n, i + 2) + +! calculate the left and right derivatives +! + dfl = f(ip1) - f(i ) + dfr = f(ip2) - f(ip1) + +! get the limited slope +! + df = limiter_minmod(0.5d+00, dfl, dfr) + ! calculate new states ! - fl(ip1) = f(ip1) - fr(i ) = f(ip1) + fl(ip1) = f(ip1) + df + fr(i ) = f(ip1) - df end if From f574dc2a330188fa6c46a48e12d7d58f8aaf2459 Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Sun, 10 May 2015 09:43:39 -0300 Subject: [PATCH 3/8] INTERPOLATIONS: Make the extrema clipping limiter independent. Signed-off-by: Grzegorz Kowal --- src/interpolations.F90 | 31 +++++++++++++++++++++++++++++-- 1 file changed, 29 insertions(+), 2 deletions(-) diff --git a/src/interpolations.F90 b/src/interpolations.F90 index 90ca9ee..bb47e86 100644 --- a/src/interpolations.F90 +++ b/src/interpolations.F90 @@ -51,6 +51,7 @@ module interpolations ! procedure(reconstruct) , pointer, save :: reconstruct_states => null() procedure(limiter_zero) , pointer, save :: limiter => null() + procedure(limiter_zero) , pointer, save :: limiter_clip => null() ! module parameters ! @@ -112,10 +113,12 @@ module interpolations ! character(len=255) :: sreconstruction = "tvd" character(len=255) :: slimiter = "mm" + character(len=255) :: climiter = "mm" character(len=255) :: positivity_fix = "off" character(len=255) :: clip_extrema = "off" character(len=255) :: name_rec = "" character(len=255) :: name_lim = "" + character(len=255) :: name_clim = "" ! !------------------------------------------------------------------------------- ! @@ -138,6 +141,7 @@ module interpolations call get_parameter_string ("limiter" , slimiter ) call get_parameter_string ("fix_positivity" , positivity_fix ) call get_parameter_string ("clip_extrema" , clip_extrema ) + call get_parameter_string ("extrema_limiter", climiter ) call get_parameter_integer("nghosts" , ng ) call get_parameter_real ("eps" , eps ) call get_parameter_real ("limo3_rad" , rad ) @@ -230,6 +234,26 @@ module interpolations limiter => limiter_zero end select +! select the clipping limiter +! + select case(trim(climiter)) + case ("mm", "minmod") + name_clim = "minmod" + limiter_clip => limiter_minmod + case ("mc", "monotonized_central") + name_clim = "monotonized central" + limiter_clip => limiter_monotonized_central + case ("sb", "superbee") + name_clim = "superbee" + limiter_clip => limiter_superbee + case ("vl", "vanleer") + name_clim = "van Leer" + limiter_clip => limiter_vanleer + case default + name_clim = "zero derivative" + limiter_clip => limiter_zero + end select + ! check additional reconstruction limiting ! select case(trim(positivity_fix)) @@ -253,6 +277,9 @@ module interpolations write (*,"(4x,a15,8x,'=',1x,a)") "limiter ", trim(name_lim) write (*,"(4x,a15,8x,'=',1x,a)") "fix positivity ", trim(positivity_fix) write (*,"(4x,a15,8x,'=',1x,a)") "clip extrema ", trim(clip_extrema) + if (clip) then + write (*,"(4x,a15,8x,'=',1x,a)") "extrema limiter", trim(name_clim) + end if end if @@ -2766,7 +2793,7 @@ module interpolations ! get the limited slope ! - df = limiter_minmod(0.5d+00, dfl, dfr) + df = limiter_clip(0.5d+00, dfl, dfr) ! calculate new states ! @@ -2790,7 +2817,7 @@ module interpolations ! get the limited slope ! - df = limiter_minmod(0.5d+00, dfl, dfr) + df = limiter_clip(0.5d+00, dfl, dfr) ! calculate new states ! From 408760bf15790eb17120ad666b5ca23ad22652fa Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Sun, 10 May 2015 10:05:22 -0300 Subject: [PATCH 4/8] INTERPOLATIONS: Make the prolongation limiter independent. Signed-off-by: Grzegorz Kowal --- src/boundaries.F90 | 24 ++++++++++---------- src/interpolations.F90 | 51 +++++++++++++++++++++++++++++++----------- src/mesh.F90 | 8 +++---- 3 files changed, 54 insertions(+), 29 deletions(-) diff --git a/src/boundaries.F90 b/src/boundaries.F90 index f1f4503..a6bdf83 100644 --- a/src/boundaries.F90 +++ b/src/boundaries.F90 @@ -5426,7 +5426,7 @@ module boundaries use coordinates , only : im , jm , km use coordinates , only : faces_dp use equations , only : nv - use interpolations , only : limiter + use interpolations , only : limiter_prol ! local variables are not implicit by default ! @@ -5489,15 +5489,15 @@ module boundaries ! dql = qn(p,i ,j,k) - qn(p,im1,j,k) dqr = qn(p,ip1,j,k) - qn(p,i ,j,k) - dqx = limiter(0.25d+00, dql, dqr) + dqx = limiter_prol(0.25d+00, dql, dqr) dql = qn(p,i,j ,k) - qn(p,i,jm1,k) dqr = qn(p,i,jp1,k) - qn(p,i,j ,k) - dqy = limiter(0.25d+00, dql, dqr) + dqy = limiter_prol(0.25d+00, dql, dqr) dql = qn(p,i,j,k ) - qn(p,i,j,km1) dqr = qn(p,i,j,kp1) - qn(p,i,j,k ) - dqz = limiter(0.25d+00, dql, dqr) + dqz = limiter_prol(0.25d+00, dql, dqr) ! calculate the derivative terms ! @@ -5713,7 +5713,7 @@ module boundaries use coordinates , only : im , jm , km use coordinates , only : edges_dp use equations , only : nv - use interpolations , only : limiter + use interpolations , only : limiter_prol ! local variables are not implicit by default ! @@ -5789,16 +5789,16 @@ module boundaries ! dql = qn(p,i ,j,k) - qn(p,im1,j,k) dqr = qn(p,ip1,j,k) - qn(p,i ,j,k) - dqx = limiter(0.25d+00, dql, dqr) + dqx = limiter_prol(0.25d+00, dql, dqr) dql = qn(p,i,j ,k) - qn(p,i,jm1,k) dqr = qn(p,i,jp1,k) - qn(p,i,j ,k) - dqy = limiter(0.25d+00, dql, dqr) + dqy = limiter_prol(0.25d+00, dql, dqr) #if NDIMS == 3 dql = qn(p,i,j,k ) - qn(p,i,j,km1) dqr = qn(p,i,j,kp1) - qn(p,i,j,k ) - dqz = limiter(0.25d+00, dql, dqr) + dqz = limiter_prol(0.25d+00, dql, dqr) #endif /* NDIMS == 3 */ #if NDIMS == 2 @@ -5969,7 +5969,7 @@ module boundaries use coordinates , only : im , jm , km use coordinates , only : corners_dp use equations , only : nv - use interpolations , only : limiter + use interpolations , only : limiter_prol ! local variables are not implicit by default ! @@ -6050,16 +6050,16 @@ module boundaries ! dql = qn(p,i ,j,k) - qn(p,im1,j,k) dqr = qn(p,ip1,j,k) - qn(p,i ,j,k) - dqx = limiter(0.25d+00, dql, dqr) + dqx = limiter_prol(0.25d+00, dql, dqr) dql = qn(p,i,j ,k) - qn(p,i,jm1,k) dqr = qn(p,i,jp1,k) - qn(p,i,j ,k) - dqy = limiter(0.25d+00, dql, dqr) + dqy = limiter_prol(0.25d+00, dql, dqr) #if NDIMS == 3 dql = qn(p,i,j,k ) - qn(p,i,j,km1) dqr = qn(p,i,j,kp1) - qn(p,i,j,k ) - dqz = limiter(0.25d+00, dql, dqr) + dqz = limiter_prol(0.25d+00, dql, dqr) #endif /* NDIMS == 3 */ #if NDIMS == 2 diff --git a/src/interpolations.F90 b/src/interpolations.F90 index bb47e86..d4dc4cd 100644 --- a/src/interpolations.F90 +++ b/src/interpolations.F90 @@ -51,6 +51,7 @@ module interpolations ! procedure(reconstruct) , pointer, save :: reconstruct_states => null() procedure(limiter_zero) , pointer, save :: limiter => null() + procedure(limiter_zero) , pointer, save :: limiter_prol => null() procedure(limiter_zero) , pointer, save :: limiter_clip => null() ! module parameters @@ -74,7 +75,7 @@ module interpolations ! declare public subroutines ! public :: initialize_interpolations, finalize_interpolations - public :: reconstruct, limiter + public :: reconstruct, limiter_prol public :: fix_positivity !- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - @@ -113,11 +114,13 @@ module interpolations ! character(len=255) :: sreconstruction = "tvd" character(len=255) :: slimiter = "mm" + character(len=255) :: plimiter = "mm" character(len=255) :: climiter = "mm" character(len=255) :: positivity_fix = "off" character(len=255) :: clip_extrema = "off" character(len=255) :: name_rec = "" character(len=255) :: name_lim = "" + character(len=255) :: name_plim = "" character(len=255) :: name_clim = "" ! !------------------------------------------------------------------------------- @@ -137,14 +140,15 @@ module interpolations ! obtain the user defined interpolation methods and coefficients ! - call get_parameter_string ("reconstruction" , sreconstruction) - call get_parameter_string ("limiter" , slimiter ) - call get_parameter_string ("fix_positivity" , positivity_fix ) - call get_parameter_string ("clip_extrema" , clip_extrema ) - call get_parameter_string ("extrema_limiter", climiter ) - call get_parameter_integer("nghosts" , ng ) - call get_parameter_real ("eps" , eps ) - call get_parameter_real ("limo3_rad" , rad ) + call get_parameter_string ("reconstruction" , sreconstruction) + call get_parameter_string ("limiter" , slimiter ) + call get_parameter_string ("fix_positivity" , positivity_fix ) + call get_parameter_string ("clip_extrema" , clip_extrema ) + call get_parameter_string ("extrema_limiter" , climiter ) + call get_parameter_string ("prolongation_limiter", plimiter ) + call get_parameter_integer("nghosts" , ng ) + call get_parameter_real ("eps" , eps ) + call get_parameter_real ("limo3_rad" , rad ) ! select the reconstruction method ! @@ -234,6 +238,26 @@ module interpolations limiter => limiter_zero end select +! select the prolongation limiter +! + select case(trim(plimiter)) + case ("mm", "minmod") + name_plim = "minmod" + limiter_prol => limiter_minmod + case ("mc", "monotonized_central") + name_plim = "monotonized central" + limiter_prol => limiter_monotonized_central + case ("sb", "superbee") + name_plim = "superbee" + limiter_prol => limiter_superbee + case ("vl", "vanleer") + name_plim = "van Leer" + limiter_prol => limiter_vanleer + case default + name_plim = "zero derivative" + limiter_prol => limiter_zero + end select + ! select the clipping limiter ! select case(trim(climiter)) @@ -273,10 +297,11 @@ module interpolations ! if (verbose) then - write (*,"(4x,a15,8x,'=',1x,a)") "reconstruction ", trim(name_rec) - write (*,"(4x,a15,8x,'=',1x,a)") "limiter ", trim(name_lim) - write (*,"(4x,a15,8x,'=',1x,a)") "fix positivity ", trim(positivity_fix) - write (*,"(4x,a15,8x,'=',1x,a)") "clip extrema ", trim(clip_extrema) + write (*,"(4x,a14, 9x,'=',1x,a)") "reconstruction" , trim(name_rec) + write (*,"(4x, a7,16x,'=',1x,a)") "limiter" , trim(name_lim) + write (*,"(4x,a20, 3x,'=',1x,a)") "prolongation limiter", trim(name_plim) + write (*,"(4x,a14, 9x,'=',1x,a)") "fix positivity" , trim(positivity_fix) + write (*,"(4x,a12,11x,'=',1x,a)") "clip extrema" , trim(clip_extrema) if (clip) then write (*,"(4x,a15,8x,'=',1x,a)") "extrema limiter", trim(name_clim) end if diff --git a/src/mesh.F90 b/src/mesh.F90 index 6f78604..ed58cdf 100644 --- a/src/mesh.F90 +++ b/src/mesh.F90 @@ -1055,7 +1055,7 @@ module mesh use coordinates , only : ng, nh, in, jn, kn, im, jm, km use coordinates , only : ib, ie, jb, je, kb, ke use equations , only : nv - use interpolations , only : limiter + use interpolations , only : limiter_prol ! local variables are not implicit by default ! @@ -1141,16 +1141,16 @@ module mesh dul = pdata%u(p,i ,j,k) - pdata%u(p,i-1,j,k) dur = pdata%u(p,i+1,j,k) - pdata%u(p,i ,j,k) - dux = limiter(0.25d+00, dul, dur) + dux = limiter_prol(0.25d+00, dul, dur) dul = pdata%u(p,i,j ,k) - pdata%u(p,i,j-1,k) dur = pdata%u(p,i,j+1,k) - pdata%u(p,i,j ,k) - duy = limiter(0.25d+00, dul, dur) + duy = limiter_prol(0.25d+00, dul, dur) #if NDIMS == 3 dul = pdata%u(p,i,j,k ) - pdata%u(p,i,j,k-1) dur = pdata%u(p,i,j,k+1) - pdata%u(p,i,j,k ) - duz = limiter(0.25d+00, dul, dur) + duz = limiter_prol(0.25d+00, dul, dur) #endif /* NDIMS == 3 */ #if NDIMS == 2 From 190f1ba6445c860bf2770039ba2beb12b1068a1d Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Sun, 10 May 2015 10:12:20 -0300 Subject: [PATCH 5/8] INTERPOLATIONS: Rename limiter() to limiter_tvd(). Signed-off-by: Grzegorz Kowal --- src/interpolations.F90 | 44 ++++++++++++++++++++++-------------------- 1 file changed, 23 insertions(+), 21 deletions(-) diff --git a/src/interpolations.F90 b/src/interpolations.F90 index d4dc4cd..131cfd6 100644 --- a/src/interpolations.F90 +++ b/src/interpolations.F90 @@ -50,7 +50,7 @@ module interpolations ! pointers to the reconstruction and limiter procedures ! procedure(reconstruct) , pointer, save :: reconstruct_states => null() - procedure(limiter_zero) , pointer, save :: limiter => null() + procedure(limiter_zero) , pointer, save :: limiter_tvd => null() procedure(limiter_zero) , pointer, save :: limiter_prol => null() procedure(limiter_zero) , pointer, save :: limiter_clip => null() @@ -113,13 +113,13 @@ module interpolations ! local variables ! character(len=255) :: sreconstruction = "tvd" - character(len=255) :: slimiter = "mm" + character(len=255) :: tlimiter = "mm" character(len=255) :: plimiter = "mm" character(len=255) :: climiter = "mm" character(len=255) :: positivity_fix = "off" character(len=255) :: clip_extrema = "off" character(len=255) :: name_rec = "" - character(len=255) :: name_lim = "" + character(len=255) :: name_tlim = "" character(len=255) :: name_plim = "" character(len=255) :: name_clim = "" ! @@ -141,7 +141,7 @@ module interpolations ! obtain the user defined interpolation methods and coefficients ! call get_parameter_string ("reconstruction" , sreconstruction) - call get_parameter_string ("limiter" , slimiter ) + call get_parameter_string ("limiter" , tlimiter ) call get_parameter_string ("fix_positivity" , positivity_fix ) call get_parameter_string ("clip_extrema" , clip_extrema ) call get_parameter_string ("extrema_limiter" , climiter ) @@ -215,27 +215,27 @@ module interpolations end if end select -! select the limiter +! select the TVD limiter ! - select case(trim(slimiter)) + select case(trim(tlimiter)) case ("mm", "minmod") - name_lim = "minmod" - limiter => limiter_minmod + name_tlim = "minmod" + limiter_tvd => limiter_minmod case ("mc", "monotonized_central") - name_lim = "monotonized central" - limiter => limiter_monotonized_central + name_tlim = "monotonized central" + limiter_tvd => limiter_monotonized_central case ("sb", "superbee") - name_lim = "superbee" - limiter => limiter_superbee + name_tlim = "superbee" + limiter_tvd => limiter_superbee case ("vl", "vanleer") - name_lim = "van Leer" - limiter => limiter_vanleer + name_tlim = "van Leer" + limiter_tvd => limiter_vanleer case ("va", "vanalbada") - name_lim = "van Albada" - limiter => limiter_vanalbada + name_tlim = "van Albada" + limiter_tvd => limiter_vanalbada case default - name_lim = "zero derivative" - limiter => limiter_zero + name_tlim = "zero derivative" + limiter_tvd => limiter_zero end select ! select the prolongation limiter @@ -298,7 +298,7 @@ module interpolations if (verbose) then write (*,"(4x,a14, 9x,'=',1x,a)") "reconstruction" , trim(name_rec) - write (*,"(4x, a7,16x,'=',1x,a)") "limiter" , trim(name_lim) + write (*,"(4x,a11,12x,'=',1x,a)") "TVD limiter" , trim(name_tlim) write (*,"(4x,a20, 3x,'=',1x,a)") "prolongation limiter", trim(name_plim) write (*,"(4x,a14, 9x,'=',1x,a)") "fix positivity" , trim(positivity_fix) write (*,"(4x,a12,11x,'=',1x,a)") "clip extrema" , trim(clip_extrema) @@ -350,7 +350,9 @@ module interpolations ! release the procedure pointers ! nullify(reconstruct_states) - nullify(limiter) + nullify(limiter_tvd) + nullify(limiter_prol) + nullify(limiter_clip) #ifdef PROFILE ! stop accounting time for module initialization/finalization @@ -471,7 +473,7 @@ module interpolations ! obtain the TVD limited derivative ! - df = limiter(0.5d+00, dfl, dfr) + df = limiter_tvd(0.5d+00, dfl, dfr) ! update the left and right-side interpolation states ! From b1f397e7f85264ef4a62c1c697ca6cfb1db38588 Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Mon, 11 May 2015 07:03:24 -0300 Subject: [PATCH 6/8] SHAPES: Calculate conserved variables for jet only once. Signed-off-by: Grzegorz Kowal --- src/shapes.F90 | 13 ++++++------- 1 file changed, 6 insertions(+), 7 deletions(-) diff --git a/src/shapes.F90 b/src/shapes.F90 index 4c47fdc..8480b20 100644 --- a/src/shapes.F90 +++ b/src/shapes.F90 @@ -601,7 +601,7 @@ module shapes ! local arrays ! real(kind=8), dimension(nv,im) :: q, u - real(kind=8), dimension(nv) :: qj + real(kind=8), dimension(nv) :: qj, uj real(kind=8), dimension(im) :: x real(kind=8), dimension(jm) :: y real(kind=8), dimension(km) :: z @@ -650,6 +650,7 @@ module shapes qj(ibz) = bjet qj(ibp) = 0.0d+00 end if ! ibx > 0 + call prim2cons(1, qj(1:nv), uj(1:nv)) ! prepare block coordinates ! @@ -674,6 +675,7 @@ module shapes ! copy the primitive variable vector ! q(1:nv,1:im) = pdata%q(1:nv,1:im,j,k) + u(1:nv,1:im) = pdata%u(1:nv,1:im,j,k) ! calculate radius ! @@ -683,22 +685,19 @@ module shapes do i = 1, im if (x(i) <= max(dx, ljet)) then q(1:nv,i) = qj(1:nv) + u(1:nv,i) = uj(1:nv) end if end do ! i = 1, im end if ! R < Rjet -! convert the primitive variables to conservative ones +! copy the primitive variables to the current block ! - call prim2cons(im, q(1:nv,1:im), u(1:nv,1:im)) + pdata%q(1:nv,1:im,j,k) = q(1:nv,1:im) ! copy the conserved variables to the current block ! pdata%u(1:nv,1:im,j,k) = u(1:nv,1:im) -! copy the primitive variables to the current block -! - pdata%q(1:nv,1:im,j,k) = q(1:nv,1:im) - end do ! j = 1, jm end do ! k = 1, km From 2e4dfbdc09cd48d2f9c92e0ea5e5e8a819af07ff Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Mon, 11 May 2015 07:33:29 -0300 Subject: [PATCH 7/8] EVOLUTION: Move boundary_variables() to update_variables(). Signed-off-by: Grzegorz Kowal --- src/evolution.F90 | 98 +++-------------------------------------------- 1 file changed, 6 insertions(+), 92 deletions(-) diff --git a/src/evolution.F90 b/src/evolution.F90 index e93af5f..b305048 100644 --- a/src/evolution.F90 +++ b/src/evolution.F90 @@ -287,7 +287,6 @@ module evolution ! include external procedures ! use blocks , only : set_blocks_update - use boundaries , only : boundary_variables use mesh , only : update_mesh ! include external variables @@ -334,10 +333,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef DEBUG ! check variables for NaNs ! @@ -516,7 +511,6 @@ module evolution ! include external procedures ! - use boundaries , only : boundary_variables use sources , only : update_sources ! include external variables @@ -589,10 +583,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef PROFILE ! stop accounting time for one step update ! @@ -623,7 +613,6 @@ module evolution ! include external procedures ! - use boundaries , only : boundary_variables use sources , only : update_sources ! include external variables @@ -691,10 +680,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - ! update fluxes from the intermediate stage ! call update_fluxes() @@ -740,10 +725,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef PROFILE ! stop accounting time for one step update ! @@ -777,7 +758,6 @@ module evolution ! include external procedures ! - use boundaries , only : boundary_variables use sources , only : update_sources ! include external variables @@ -901,10 +881,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - end do ! n = 1, stages - 1 != the final step: U(n+1) = 1/m U(0) + (m-1)/m [1 + dt/(m-1) L] U(m-1) @@ -954,10 +930,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef PROFILE ! stop accounting time for one step update ! @@ -989,7 +961,6 @@ module evolution ! include external procedures ! - use boundaries , only : boundary_variables use sources , only : update_sources ! include external variables @@ -1072,10 +1043,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - !! 2nd substep of integration !! ! prepare the fractional time step @@ -1118,10 +1085,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - !! 3rd substep of integration !! ! prepare the fractional time step @@ -1173,10 +1136,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef PROFILE ! stop accounting time for one step update ! @@ -1209,7 +1168,6 @@ module evolution ! include external procedures ! - use boundaries , only : boundary_variables use sources , only : update_sources ! include external variables @@ -1292,10 +1250,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != 2nd step: U(2) = U(1) + 1/2 dt F[U(1)] ! ! update fluxes @@ -1333,10 +1287,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != 3rd step: U(3) = 2/3 U(n) + 1/3 U(2) + 1/6 dt F[U(2)] ! ! calculate the fractional time step @@ -1379,10 +1329,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != the final step: U(n+1) = U(3) + 1/2 dt F[U(3)] ! ! calculate fractional time step @@ -1433,10 +1379,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef PROFILE ! stop accounting time for one step update ! @@ -1469,7 +1411,6 @@ module evolution ! include external procedures ! - use boundaries , only : boundary_variables use sources , only : update_sources ! include external variables @@ -1560,10 +1501,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != 2nd step: U(2) = U(1) + b1 dt F[U(1)] ! ! update fluxes @@ -1601,10 +1538,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != 3rd step: U(3) = a31 U(n) + a33 U(2) + b3 dt F[U(2)] ! ! calculate the fractional time step @@ -1647,10 +1580,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != 4th step: U(4) = a41 U(n) + a44 U(3) + b4 dt F[U(3)] ! ! calculate the fractional time step @@ -1697,10 +1626,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - != the final step: U(n+1) = a53 U(2) + a55 U(4) + b5 dt F[U(4)] ! ! calculate the fractional time step @@ -1748,10 +1673,6 @@ module evolution ! call update_variables() -! update boundaries -! - call boundary_variables() - #ifdef PROFILE ! stop accounting time for one step update ! @@ -1961,6 +1882,7 @@ module evolution ! include external procedures ! + use boundaries , only : boundary_variables use equations , only : update_primitive_variables use shapes , only : update_shapes @@ -1986,32 +1908,24 @@ module evolution call start_timer(imv) #endif /* PROFILE */ -! associate the pointer with the first block on the data block list +! update primitive variables and shapes if necessary ! pdata => list_data - -! iterate over all data blocks -! do while (associated(pdata)) - -! associate pmeta with the corresponding meta block -! pmeta => pdata%meta -! convert conserved variables to primitive ones for the current block and -! update shapes if necessary -! if (pmeta%update) then call update_primitive_variables(pdata%u, pdata%q) call update_shapes(pdata) end if -! assign pointer to the next block -! pdata => pdata%next - end do +! update boundaries +! + call boundary_variables() + #ifdef PROFILE ! stop accounting time for variable update ! From 8bc24eeb80a7c2f2156c60ae5fadbd82add62674 Mon Sep 17 00:00:00 2001 From: Grzegorz Kowal Date: Mon, 11 May 2015 10:38:09 -0300 Subject: [PATCH 8/8] SCHEMES: In HD and MHD we reconstruct states using primitive variables. Signed-off-by: Grzegorz Kowal --- src/schemes.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/schemes.F90 b/src/schemes.F90 index a0eb0c8..d6ed6ad 100644 --- a/src/schemes.F90 +++ b/src/schemes.F90 @@ -114,7 +114,7 @@ module schemes character(len=64) :: solver = "HLL" character(len=64) :: statev = "primitive" character(len=255) :: name_sol = "" - character(len=255) :: name_sts = "" + character(len=255) :: name_sts = "primitive" ! !------------------------------------------------------------------------------- !