diff --git a/Source/radiation/Make.package b/Source/radiation/Make.package index f564af29cf..fb5c5990d8 100644 --- a/Source/radiation/Make.package +++ b/Source/radiation/Make.package @@ -42,7 +42,6 @@ ca_F90EXE_sources += rad_params.F90 ca_F90EXE_sources += Rad_nd.F90 CEXE_headers += fluxlimiter.H CEXE_headers += RadHydro.H -ca_F90EXE_sources += fluxlimiter.F90 ca_F90EXE_sources += filter.F90 CEXE_sources += RadDerive.cpp CEXE_headers += RadDerive.H diff --git a/Source/radiation/RAD_F.H b/Source/radiation/RAD_F.H index 2306d83fa6..55e0154615 100644 --- a/Source/radiation/RAD_F.H +++ b/Source/radiation/RAD_F.H @@ -65,16 +65,6 @@ extern "C" { } #endif -// -#ifdef __cplusplus -extern "C" { -#endif -void ca_initfluxlimiter - (const int* limiter, const int* closure); -#ifdef __cplusplus -} -#endif - extern "C" { void ca_rad_source diff --git a/Source/radiation/Radiation.cpp b/Source/radiation/Radiation.cpp index ace2a6969b..d2eb374b1b 100644 --- a/Source/radiation/Radiation.cpp +++ b/Source/radiation/Radiation.cpp @@ -333,8 +333,6 @@ Radiation::Radiation(Amr* Parent, Castro* castro, int restart) amrex::Abort("MGFLDSolver does not support limiter = 1"); } - ca_initfluxlimiter(&radiation::limiter, &radiation::closure); - inner_update_limiter = 0; pp.query("inner_update_limiter", inner_update_limiter); diff --git a/Source/radiation/fluxlimiter.F90 b/Source/radiation/fluxlimiter.F90 deleted file mode 100644 index 1c29118b36..0000000000 --- a/Source/radiation/fluxlimiter.F90 +++ /dev/null @@ -1,140 +0,0 @@ -module fluxlimiter_module - - use amrex_fort_module, only: rt => amrex_real - - implicit none - - integer, allocatable, save :: limiter, closure - - real(rt), parameter :: OneThird = 1.e0_rt/3.e0_rt, TwoThirds=2.e0_rt/3.e0_rt, & - OneSixth=1.e0_rt/6.e0_rt, TwoNinths=2.e0_rt/9.e0_rt, & - FiveNinths=5.e0_rt/9.e0_rt - -contains - - subroutine init_fluxlimiter_module(limiter_in, closure_in) - - use amrex_fort_module, only: rt => amrex_real - - implicit none - - integer, intent(in) :: limiter_in, closure_in - - allocate(limiter) - allocate(closure) - - limiter = limiter_in - closure = closure_in - - end subroutine init_fluxlimiter_module - - function FLDlambda(r, limiter_in) - - use amrex_fort_module, only: rt => amrex_real - - integer, intent(in), optional :: limiter_in - real(rt), intent(in) :: r - - real(rt) :: FLDlambda - integer :: limiter_local - - if (present(limiter_in))then - limiter_local = limiter_in - else - limiter_local = limiter - end if - - if (limiter_local .eq. 0) then ! no limiter - FLDlambda = OneThird - else if (limiter_local < 10) then ! approximate LP, [123] - FLDlambda = (2.e0_rt + r) / (6.e0_rt + r * (3.e0_rt + r)) - else if (limiter_local < 20) then ! Bruenn, 1[123] - FLDlambda = 1.e0_rt / (3.e0_rt + r) - else if (limiter_local < 30) then ! Larsen's square root, 2[123] - FLDlambda = 1.e0_rt / sqrt(9.e0_rt + r**2) - else if (limiter_local < 40) then ! Minerbo, 3[123] - if (r .lt. 1.5e0_rt) then - FLDlambda = 2.e0_rt/(3.e0_rt + sqrt(9.e0_rt+12.e0_rt*r**2)) - else - FLDlambda = 1.e0_rt/(1.e0_rt+r+sqrt(1.e0_rt+2.e0_rt*r)) - end if - else - print *, "Unknown limiter ", limiter_local - stop - endif - - end function FLDlambda - - function Edd_factor(lambda) result(f) - - use amrex_fort_module, only: rt => amrex_real - - implicit none - - real(rt), intent(in) :: lambda - real(rt) :: f - - if (closure .eq. 0) then - f = lambda - else if (closure .eq. 1) then - f = OneThird - else if (closure .eq. 2) then - f = 1.e0_rt - 2.e0_rt * lambda - else if (closure .eq. 3) then ! lambda + (lambda*R)**2 - if (limiter .eq. 0) then ! no limiter - f = OneThird - else if (limiter < 10) then ! approximate LP, [123] - f = lambda + 0.25e0_rt*(max(0.0_rt, (1.e0_rt-3.e0_rt*lambda)) + & - sqrt(max(0.0_rt, (1.e0_rt-3.e0_rt*lambda))*(1.e0_rt+5.e0_rt*lambda)))**2 - else if (limiter < 20) then ! Bruenn, 1[123] - f = 1.0e0_rt - 5.e0_rt*lambda + 9.e0_rt*lambda**2 - else if (limiter < 30) then ! Larsen's square root, 2[123] - f = 1.0e0_rt + lambda - 9.e0_rt*lambda**2 - else if (limiter < 40) then ! Minerbo - if (lambda .gt. TwoNinths) then - f = OneThird - else - f = 1.e0_rt + 3.e0_rt*lambda - 2.e0_rt*sqrt(2.e0_rt*lambda) - end if - else - print *, "Unknown limiter ", limiter - stop - endif - else if (closure .eq. 4) then ! 1/3 + 2/3*(lambda*R)**2 - if (limiter .eq. 0) then ! no limiter - f = OneThird - else if (limiter < 10) then ! approximate LP, [123] - f = OneThird + OneSixth*(max(0.0_rt, (1.e0_rt-3.e0_rt*lambda)) + & - sqrt(max(0.0_rt, (1.e0_rt-3.e0_rt*lambda))*(1.e0_rt+5.e0_rt*lambda)))**2 - else if (limiter < 20) then ! Bruenn, 1[123] - f = OneThird + TwoThirds*(1.0e0_rt - 6.e0_rt*lambda + 9.e0_rt*lambda**2) - else if (limiter < 30) then ! Larsen's square root, 2[123] - f = OneThird + TwoThirds*(1.0e0_rt - 9.e0_rt*lambda**2) - else if (limiter < 40) then ! Minerbo - if (lambda .gt. TwoNinths) then - f = FiveNinths - TwoThirds*lambda - else - f = OneThird + TwoThirds*(1.e0_rt + 2.e0_rt*lambda - 2.e0_rt*sqrt(2.e0_rt*lambda)) - end if - else - print *, "Unknown limiter ", limiter - stop - endif - end if - - end function Edd_factor - -end module fluxlimiter_module - -subroutine ca_initfluxlimiter(limiter, closure) bind(C, name="ca_initfluxlimiter") - - use fluxlimiter_module, only: init_fluxlimiter_module - use amrex_fort_module, only: rt => amrex_real - - implicit none - - integer, intent(in) :: limiter, closure - - call init_fluxlimiter_module(limiter, closure) - -end subroutine ca_initfluxlimiter