33 double precision :: sts_dtpar=0.8d0
36 double precision,
parameter :: nu_sts = 0.5d0
38 integer :: sts_ncycles=1000
39 integer :: sts_method = 1
45 logical :: fix_conserve_at_step = .true.
46 logical :: sts_initialized = .false.
51 subroutine subr1(ixI^L,ixO^L,w,x,wres,fix_conserve_at_step,my_dt,igrid,nflux)
53 integer,
intent(in) :: ixi^
l, ixo^
l, igrid, nflux
54 double precision,
intent(in) :: x(ixi^s,1:
ndim)
55 double precision,
intent(inout) :: wres(ixi^s,1:nw), w(ixi^s,1:nw)
56 double precision,
intent(in) :: my_dt
57 logical,
intent(in) :: fix_conserve_at_step
61 function subr2(w,ixG^L,ix^L,dx^D,x)
result(dtnew)
63 integer,
intent(in) :: ixg^
l, ix^
l
64 double precision,
intent(in) ::
dx^
d, x(ixg^s,1:
ndim)
65 double precision,
intent(in) :: w(ixg^s,1:nw)
66 double precision :: dtnew
70 subroutine subr_e(w, x, ixI^L, ixO^L, step)
73 integer,
intent(in) :: ixi^
l,ixo^
l
74 double precision,
intent(inout) :: w(ixi^s,1:nw)
75 double precision,
intent(in) :: x(ixi^s,1:
ndim)
76 integer,
intent(in) :: step
80 subroutine subr5(ixI^L, ixO^L, w, x)
82 integer,
intent(in) :: ixi^
l, ixo^
l
83 double precision,
intent(in) :: x(ixi^s,1:
ndim)
84 double precision,
intent(inout) :: w(ixi^s,1:nw)
91 double precision,
intent(in) :: dt
95 function subr4(dt,dtnew,dt_modified)
result(s)
96 double precision,
intent(in) :: dtnew
97 double precision,
intent(inout) :: dt
98 logical,
intent(inout) :: dt_modified
106 double precision :: dt_expl
110 integer,
dimension(-1:1^D&) :: type_send_srl_sts_1, type_recv_srl_sts_1
111 integer,
dimension(-1:1^D&) :: type_send_r_sts_1
112 integer,
dimension( 0:3^D&) :: type_recv_r_sts_1
113 integer,
dimension( 0:3^D&) :: type_recv_p_sts_1, type_send_p_sts_1
115 integer,
dimension(-1:1^D&) :: type_send_srl_sts_2, type_recv_srl_sts_2
116 integer,
dimension(-1:1^D&) :: type_send_r_sts_2
117 integer,
dimension( 0:3^D&) :: type_recv_r_sts_2
118 integer,
dimension( 0:3^D&) :: type_recv_p_sts_2, type_send_p_sts_2
126 logical :: types_initialized
127 logical :: evolve_magnetic_field
128 procedure(subr1),
pointer,
nopass :: sts_set_sources
129 procedure(subr2),
pointer,
nopass :: sts_getdt
130 procedure(subr5),
pointer,
nopass :: sts_before_first_cycle, sts_after_last_cycle
131 procedure(subr_e),
pointer,
nopass :: sts_handle_errors
132 type(sts_term),
pointer :: next
140 procedure(subr4),
pointer :: sts_get_ncycles
148 if(.not. sts_initialized)
then
151 sts_dtpar=sts_dtpar/dble(
ndim)
152 sts_initialized = .true.
153 if(sts_method .eq. 1)
then
155 sts_get_ncycles => sts_get_ncycles1
156 else if(sts_method .eq. 2)
then
158 sts_get_ncycles => sts_get_ncycles2
160 call mpistop(
"Unknown sts method")
168 if (sts_initialized)
then
176 subroutine sts_params_read(files)
178 character(len=*),
intent(in) :: files(:)
181 namelist /sts_list/ sts_dtpar,sts_ncycles,sts_method,
sourcetype_sts
183 do n = 1,
size(files)
184 open(
unitpar, file=trim(files(n)), status=
"old")
185 read(
unitpar, sts_list,
end=111)
189 end subroutine sts_params_read
198 subroutine add_sts_method(sts_getdt, sts_set_sources, startVar, nflux, startwbc, nwbc, evolve_B)
202 integer,
intent(in) :: startvar, nflux, startwbc, nwbc
203 logical,
intent(in) :: evolve_b
207 subroutine sts_set_sources(ixI^L,ixO^L,w,x,wres,fix_conserve_at_step,my_dt,igrid,nflux)
210 integer,
intent(in) :: ixi^
l, ixo^
l, igrid, nflux
211 double precision,
intent(in) :: x(ixi^s,1:
ndim)
212 double precision,
intent(inout) :: wres(ixi^s,1:nw), w(ixi^s,1:nw)
213 double precision,
intent(in) :: my_dt
214 logical,
intent(in) :: fix_conserve_at_step
215 end subroutine sts_set_sources
217 function sts_getdt(w,ixG^L,ix^L,dx^D,x)
result(dtnew)
219 integer,
intent(in) :: ixg^
l, ix^
l
220 double precision,
intent(in) ::
dx^
d, x(ixg^s,1:
ndim)
221 double precision,
intent(in) :: w(ixg^s,1:nw)
222 double precision :: dtnew
223 end function sts_getdt
227 type(sts_term),
pointer :: temp
230 temp%sts_getdt => sts_getdt
231 temp%sts_set_sources => sts_set_sources
232 temp%sts_before_first_cycle => null()
233 temp%sts_after_last_cycle => null()
234 temp%sts_handle_errors => null()
235 temp%startVar = startvar
236 temp%endVar= startvar+nflux-1
238 temp%startwbc = startwbc
240 temp%types_initialized = .false.
241 temp%evolve_magnetic_field=evolve_b
253 subroutine sts_before_first_cycle(ixI^L, ixO^L, w, x)
255 integer,
intent(in) :: ixi^
l, ixo^
l
256 double precision,
intent(in) :: x(ixi^s,1:
ndim)
257 double precision,
intent(inout) :: w(ixi^s,1:nw)
258 end subroutine sts_before_first_cycle
260 subroutine sts_after_last_cycle(ixI^L, ixO^L, w, x)
262 integer,
intent(in) :: ixi^
l, ixo^
l
263 double precision,
intent(in) :: x(ixi^s,1:
ndim)
264 double precision,
intent(inout) :: w(ixi^s,1:nw)
265 end subroutine sts_after_last_cycle
278 subroutine sts_error_handling(w, x, ixI^L, ixO^L, step)
281 integer,
intent(in) :: ixi^
l,ixo^
l
282 double precision,
intent(inout) :: w(ixi^s,1:nw)
283 double precision,
intent(in) :: x(ixi^s,1:
ndim)
284 integer,
intent(in) :: step
285 end subroutine sts_error_handling
292 function sts_get_ncycles1(dt,dtnew,dt_modified)
result(is)
293 double precision,
intent(in) :: dtnew
294 double precision,
intent(inout) :: dt
295 logical,
intent(inout) :: dt_modified
298 double precision :: ss
301 ss = dtnew*((2.d0*sts_ncycles+1)**2-9.d0)/16.d0
309 if(ss .le. 1.d0)
then
312 is=ceiling((dsqrt(9.d0+16.d0*ss)-1.d0)*0.5d0)
317 end function sts_get_ncycles1
320 function sts_get_ncycles2(dt,dtnew,dt_modified)
result(is)
321 double precision,
intent(in) :: dtnew
322 double precision,
intent(inout) :: dt
323 logical,
intent(inout) :: dt_modified
326 double precision :: ss,rr
331 ncycles = sts_ncycles
333 ss=sum_chev(nu_sts,ncycles,rr)
342 end function sts_get_ncycles2
348 double precision,
intent(inout) :: my_dt
349 double precision :: my_dt1
350 logical :: dt_modified, dt_modified1, dt_modified2
352 double precision :: dtnew,dtmin_mype
353 double precision ::
dx^
d, ss
354 integer:: iigrid, igrid, ncycles
355 type(sts_term),
pointer :: temp,oldtemp
358 dt_modified = .false.
359 do while(
associated(temp))
360 dt_modified2 = .false.
363 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
370 dtmin_mype=min(dtmin_mype, sts_dtpar * temp%sts_getdt(ps(igrid)%w,ixg^
ll,
ixm^
ll,
dx^
d,ps(igrid)%x))
373 call mpi_allreduce(dtmin_mype,dtnew,1,mpi_double_precision,mpi_min,
icomm,
ierrmpi)
374 temp%s = sts_get_ncycles(my_dt,dtnew,dt_modified2)
382 if(dt_modified2)
then
387 dt_modified1 = .false.
388 do while(.not.
associated(oldtemp,temp))
389 oldtemp%s = sts_get_ncycles(my_dt1,oldtemp%dt_expl,dt_modified1)
391 if(dt_modified1)
call mpistop(
"sts dt modified twice")
392 oldtemp=>oldtemp%next
401 pure FUNCTION chev(j,nu,N)
404 double precision,
INTENT(IN) :: nu
405 INTEGER,
INTENT(IN) :: j, n
406 double precision :: chev
408 chev = 1d0 / ((-1d0 + nu)*cos(((2d0*j - 1d0) / n)* (
dpi/2d0)) + 1d0 + nu)
412 FUNCTION sum_chev(nu,N,limMax)
413 double precision,
intent(in) :: nu,limmax
414 integer,
intent(inout) :: n
415 double precision :: sum_chev, tmp
421 do while (j < n .and. sum_chev < limmax)
422 sum_chev = sum_chev + chev(j,nu,n)
426 END FUNCTION sum_chev
441 subroutine sts_add_source2(my_dt)
448 double precision,
intent(in) :: my_dt
449 double precision,
allocatable :: bj(:)
450 double precision :: sumbj,dtj
452 integer:: iigrid, igrid, j, ixc^
l
453 logical :: stagger_flag=.false., prolong_flag=.false., coarsen_flag=.false.
454 type(sts_term),
pointer :: temp
462 do while(
associated(temp))
464 if(.not.temp%evolve_magnetic_field)
then
475 if(
associated(temp%sts_before_first_cycle))
then
480 do iigrid=1,igridstail; igrid=igrids(iigrid);
483 call temp%sts_before_first_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
487 allocate(bj(1:temp%s))
489 bj(j) = chev(j,nu_sts,sts_ncycles)
499 if(.not. temp%types_initialized)
then
501 if(temp%nflux>temp%nwbc)
then
517 temp%types_initialized = .true.
522 if(j .eq. temp%s .and. (sumbj + bj(j)) * temp%dt_expl > my_dt)
then
523 dtj = my_dt - sumbj * temp%dt_expl
525 dtj = bj(j)* temp%dt_expl
527 sumbj = sumbj + bj(j)
530 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
533 call temp%sts_set_sources(ixg^
ll,
ixm^
ll,ps(igrid)%w,ps(igrid)%x,ps1(igrid)%w,fix_conserve_at_step,dtj,igrid,temp%nflux)
534 if(temp%nflux>
ndir)
then
535 ps(igrid)%w(
ixm^t,temp%startVar)=ps(igrid)%w(
ixm^t,temp%startVar)+dtj*ps1(igrid)%w(
ixm^t,temp%startVar)
537 ps(igrid)%ws(ixc^s,1:nws)=ps(igrid)%ws(ixc^s,1:nws)+dtj*ps1(igrid)%w(ixc^s,iw_mag(1:nws))
543 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
546 call temp%sts_set_sources(ixg^
ll,
ixm^
ll,ps(igrid)%w,ps(igrid)%x,ps1(igrid)%w,fix_conserve_at_step,dtj,igrid,temp%nflux)
547 ps(igrid)%w(
ixm^t,temp%startVar:temp%endVar)=ps(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+&
548 dtj*ps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
553 if(fix_conserve_at_step)
then
561 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
567 if(
associated(temp%sts_handle_errors))
then
569 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
572 call temp%sts_handle_errors(ps(igrid)%w,ps(igrid)%x,ixg^
ll,
ixm^
ll,j)
577 if(temp%nflux>temp%nwbc.and.temp%s==j)
then
591 if(
associated(temp%sts_after_last_cycle))
then
592 do iigrid=1,igridstail; igrid=igrids(iigrid);
595 call temp%sts_after_last_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
602 if(.not.temp%evolve_magnetic_field)
then
625 do iigrid=1,igridstail; igrid=igrids(iigrid);
631 end subroutine sts_add_source2
635 subroutine sts_add_source1(my_dt)
643 double precision,
intent(in) :: my_dt
644 double precision :: dtj
645 double precision :: omega1,cmu,cmut,cnu,cnut,one_mu_nu
646 double precision,
allocatable :: bj(:)
647 integer:: iigrid, igrid, j, ixc^
l, ixgext^
l
648 logical :: evenstep, stagger_flag=.false., prolong_flag=.false., coarsen_flag=.false., total_energy_flag=.true.
649 type(sts_term),
pointer :: temp
650 type(state),
dimension(:),
pointer :: tmpps1, tmpps2
658 do while(
associated(temp))
660 if(.not.temp%evolve_magnetic_field)
then
671 if(
associated(temp%sts_before_first_cycle))
then
679 do iigrid=1,igridstail; igrid=igrids(iigrid);
682 call temp%sts_before_first_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
683 if(.not.
allocated(ps2(igrid)%w))
allocate(ps2(igrid)%w(ixg^t,1:nw))
684 if(.not.
allocated(ps3(igrid)%w))
allocate(ps3(igrid)%w(ixg^t,1:nw))
685 if(.not.
allocated(ps4(igrid)%w))
allocate(ps4(igrid)%w(ixg^t,1:nw))
686 ps1(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
687 ps2(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
692 ixgext^
l=ixg^
ll^ladd1;
694 do iigrid=1,igridstail; igrid=igrids(iigrid);
695 if(.not.
allocated(ps2(igrid)%w))
then
698 if(.not.
allocated(ps3(igrid)%w))
allocate(ps3(igrid)%w(ixg^t,1:nw))
699 if(.not.
allocated(ps4(igrid)%w))
allocate(ps4(igrid)%w(ixg^t,1:nw))
700 ps1(igrid)%w=ps(igrid)%w
701 ps2(igrid)%w=ps(igrid)%w
702 ps1(igrid)%ws=ps(igrid)%ws
703 ps2(igrid)%ws=ps(igrid)%ws
708 do iigrid=1,igridstail; igrid=igrids(iigrid);
709 if(.not.
allocated(ps2(igrid)%w))
allocate(ps2(igrid)%w(ixg^t,1:nw))
710 if(.not.
allocated(ps3(igrid)%w))
allocate(ps3(igrid)%w(ixg^t,1:nw))
711 if(.not.
allocated(ps4(igrid)%w))
allocate(ps4(igrid)%w(ixg^t,1:nw))
712 ps1(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
713 ps2(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
719 allocate(bj(0:temp%s))
723 omega1=4.d0/dble(temp%s**2+temp%s-2)
737 if(.not. temp%types_initialized)
then
739 if(temp%nflux>temp%nwbc)
then
755 temp%types_initialized = .true.
760 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
764 call temp%sts_set_sources(ixg^
ll,
ixm^
ll,ps(igrid)%w,ps(igrid)%x,ps4(igrid)%w,fix_conserve_at_step,dtj,igrid,temp%nflux)
766 ps3(igrid)%w(ixc^s,temp%startVar:temp%endVar) = my_dt * ps4(igrid)%w(ixc^s,temp%startVar:temp%endVar)
767 if(temp%nflux>
ndir)
then
768 ps1(igrid)%w(
ixm^t,temp%startVar) = ps1(igrid)%w(
ixm^t,temp%startVar) + cmut * ps3(igrid)%w(
ixm^t,temp%startVar)
770 ps1(igrid)%ws(ixc^s,1:nws) = ps1(igrid)%ws(ixc^s,1:nws) + cmut * ps3(igrid)%w(ixc^s,iw_mag(1:nws))
776 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
779 call temp%sts_set_sources(ixg^
ll,
ixm^
ll,ps(igrid)%w,ps(igrid)%x,ps4(igrid)%w,fix_conserve_at_step,dtj,igrid,temp%nflux)
781 ps3(igrid)%w(
ixm^t,temp%startVar:temp%endVar) = my_dt * ps4(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
782 ps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar) = ps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar) + &
783 cmut * ps3(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
787 if(fix_conserve_at_step)
then
795 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
802 if(
associated(temp%sts_handle_errors))
then
804 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
807 call temp%sts_handle_errors(ps1(igrid)%w,ps1(igrid)%x,ixg^
ll,
ixm^
ll,1)
811 if(temp%nflux>temp%nwbc.and.temp%s==1)
then
830 bj(j)=dble(j**2+j-2)/dble(2*j*(j+1))
831 cmu=dble(2*j-1)/dble(j)*bj(j)/bj(j-1)
833 cnu=dble(1-j)/dble(j)*bj(j)/bj(j-2)
834 cnut=(bj(j-1)-1.d0)*cmut
835 one_mu_nu=1.d0-cmu-cnu
847 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
853 call temp%sts_set_sources(ixg^
ll,
ixm^
ll,tmpps1(igrid)%w,ps(igrid)%x,ps4(igrid)%w,fix_conserve_at_step,dtj,igrid,temp%nflux)
854 if(temp%nflux>
ndir)
then
855 tmpps2(igrid)%w(
ixm^t,temp%startVar)=cmu*tmpps1(igrid)%w(
ixm^t,temp%startVar)+&
856 cnu*tmpps2(igrid)%w(
ixm^t,temp%startVar)+one_mu_nu*ps(igrid)%w(
ixm^t,temp%startVar)+&
857 dtj*ps4(igrid)%w(
ixm^t,temp%startVar)+cnut*ps3(igrid)%w(
ixm^t,temp%startVar)
859 tmpps2(igrid)%ws(ixc^s,1:nws)=cmu*tmpps1(igrid)%ws(ixc^s,1:nws)+&
860 cnu*tmpps2(igrid)%ws(ixc^s,1:nws)+one_mu_nu*ps(igrid)%ws(ixc^s,1:nws)+&
861 dtj*ps4(igrid)%w(ixc^s,iw_mag(1:nws))+cnut*ps3(igrid)%w(ixc^s,iw_mag(1:nws))
867 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
873 call temp%sts_set_sources(ixg^
ll,
ixm^
ll,tmpps1(igrid)%w,ps(igrid)%x,ps4(igrid)%w,fix_conserve_at_step,dtj,igrid,temp%nflux)
874 tmpps2(igrid)%w(
ixm^t,temp%startVar:temp%endVar)=cmu*tmpps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+&
875 cnu*tmpps2(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+one_mu_nu*ps(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+&
876 dtj*ps4(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+cnut*ps3(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
880 if(fix_conserve_at_step)
then
888 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
894 if(
associated(temp%sts_handle_errors))
then
896 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
899 call temp%sts_handle_errors(tmpps2(igrid)%w,ps(igrid)%x,ixg^
ll,
ixm^
ll,j)
903 if(temp%nflux>temp%nwbc.and.temp%s==j)
then
915 evenstep=.not.evenstep
918 if(
associated(temp%sts_after_last_cycle))
then
920 do iigrid=1,igridstail; igrid=igrids(iigrid);
923 ps(igrid)%w(ixg^t,temp%startVar:temp%endVar)=tmpps2(igrid)%w(ixg^t,temp%startVar:temp%endVar)
924 call temp%sts_after_last_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
933 do iigrid=1,igridstail; igrid=igrids(iigrid);
934 ps(igrid)%w(ixg^t,temp%startVar:temp%endVar)=tmpps2(igrid)%w(ixg^t,temp%startVar:temp%endVar)
935 ps(igrid)%ws=tmpps2(igrid)%ws
940 do iigrid=1,igridstail; igrid=igrids(iigrid);
941 ps(igrid)%w(ixg^t,temp%startVar:temp%endVar)=tmpps2(igrid)%w(ixg^t,temp%startVar:temp%endVar)
949 if(.not.temp%evolve_magnetic_field)
then
972 do iigrid=1,igridstail; igrid=igrids(iigrid);
978 end subroutine sts_add_source1
subroutine, public alloc_state(igrid, s, ixgl, ixgextl, alloc_once_for_ps)
allocate memory to physical state of igrid node
Module for physical and numeric constants.
double precision, parameter dpi
Pi.
Module for flux conservation near refinement boundaries.
subroutine, public init_comm_fix_conserve(idimlim, nwfluxin)
subroutine, public fix_edges(psuse, idimlim)
subroutine, public recvflux(idimlim)
subroutine, public sendflux(idimlim)
subroutine, public fix_conserve(psb, idimlim, nw0, nwfluxin)
Module with geometry-related routines (e.g., divergence, curl)
update ghost cells of all blocks including physical boundaries
integer, dimension( :^d &), pointer type_recv_r
integer, dimension(^nd, 0:3) l
subroutine getbc(time, qdt, psb, nwstart, nwbc)
do update ghost cells of all blocks including physical boundaries
integer, dimension(0:3^d &), target type_send_p_f
integer, dimension( :^d &), pointer type_send_p
integer, dimension( :^d &), pointer type_send_srl
subroutine create_bc_mpi_datatype(nwstart, nwbc)
integer, dimension(0:3^d &), target type_recv_r_f
integer, dimension(-1:1^d &), target type_recv_srl_f
integer, dimension(-1:1^d &), target type_send_r_f
integer, dimension(-1:1^d &), target type_send_srl_f
integer, dimension( :^d &), pointer type_send_r
integer, dimension( :^d &), pointer type_recv_p
integer, dimension(0:3^d &), target type_recv_p_f
integer, dimension( :^d &), pointer type_recv_srl
This module contains definitions of global parameters and variables and some generic functions/subrou...
type(state), pointer block
Block pointer for using one block and its previous state.
integer, parameter unitpar
file handle for IO
double precision global_time
The global simulation time.
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical stagger_grid
True for using stagger grid.
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
integer icomm
The MPI communicator.
logical coarsenprimitive
coarsen primitive variables in level-jump ghost cells
double precision, dimension(:), allocatable, parameter d
integer ndir
Number of spatial dimensions (components) for vector variables.
integer ixm
the mesh range of a physical block without ghost cells
integer ierrmpi
A global MPI error return code.
logical time_advance
do time evolving
logical prolongprimitive
prolongate primitive variables in level-jump ghost cells
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
double precision, dimension(:,:), allocatable dx
double precision, dimension(^nd) dxlevel
store unstretched cell size of current level
This module defines the procedures of a physics module. It contains function pointers for the various...
logical phys_partial_ionization
Solve partially ionized one-fluid plasma.
logical phys_total_energy
Solve total energy equation or not.
procedure(sub_update_temperature), pointer phys_update_temperature
procedure(sub_face_to_center), pointer phys_face_to_center
Module for handling problematic values in simulations, such as negative pressures.
Generic supertimestepping method 1) in amrvac.par in sts_list set the following parameters which have...
integer, public sourcetype_sts
pure logical function, public is_sts_initialized()
integer, parameter, public sourcetype_sts_prior
logical function, public set_dt_sts_ncycles(my_dt)
This sets the explicit dt and calculates the number of cycles for each of the terms implemented with ...
subroutine, public add_sts_method(sts_getdt, sts_set_sources, startvar, nflux, startwbc, nwbc, evolve_b)
subroutine which added programatically a term to be calculated using STS Params: sts_getdt function c...
type(sts_term), pointer head_sts_terms
subroutine, public set_conversion_methods_to_head(sts_before_first_cycle, sts_after_last_cycle)
Set the hooks called before the first cycle and after the last cycle in the STS update This method sh...
integer, parameter, public sourcetype_sts_split
subroutine, public set_error_handling_to_head(sts_error_handling)
Set the hook of error handling in the STS update. This method is called before updating the BC....
subroutine, public sts_init()
Initialize sts module.
integer, parameter, public sourcetype_sts_after
procedure(subr3), pointer, public sts_add_source