34 double precision :: sts_dtpar=0.8d0
37 double precision,
parameter :: nu_sts = 0.5d0
39 integer :: sts_ncycles=1000
40 integer :: sts_method = 1
46 logical :: fix_conserve_at_step = .true.
47 logical :: sts_initialized = .false.
52 subroutine subr1(ixI^L,ixO^L,w,x,wres,fix_conserve_at_step,my_dt,igrid,nflux)
54 integer,
intent(in) :: ixi^
l, ixo^
l, igrid, nflux
55 double precision,
intent(in) :: x(ixi^s,1:
ndim)
56 double precision,
intent(inout) :: wres(ixi^s,1:nw), w(ixi^s,1:nw)
57 double precision,
intent(in) :: my_dt
58 logical,
intent(in) :: fix_conserve_at_step
62 function subr2(w,ixG^L,ix^L,dx^D,x)
result(dtnew)
64 integer,
intent(in) :: ixg^
l, ix^
l
65 double precision,
intent(in) ::
dx^
d, x(ixg^s,1:
ndim)
66 double precision,
intent(in) :: w(ixg^s,1:nw)
67 double precision :: dtnew
71 subroutine subr_e(w, x, ixI^L, ixO^L, step)
74 integer,
intent(in) :: ixi^
l,ixo^
l
75 double precision,
intent(inout) :: w(ixi^s,1:nw)
76 double precision,
intent(in) :: x(ixi^s,1:
ndim)
77 integer,
intent(in) :: step
81 subroutine subr5(ixI^L, ixO^L, w, x)
83 integer,
intent(in) :: ixi^
l, ixo^
l
84 double precision,
intent(in) :: x(ixi^s,1:
ndim)
85 double precision,
intent(inout) :: w(ixi^s,1:nw)
92 double precision,
intent(in) :: dt
96 function subr4(dt,dtnew,dt_modified)
result(s)
97 double precision,
intent(in) :: dtnew
98 double precision,
intent(inout) :: dt
99 logical,
intent(inout) :: dt_modified
107 double precision :: dt_expl
111 integer,
dimension(-1:1^D&) :: type_send_srl_sts_1, type_recv_srl_sts_1
112 integer,
dimension(-1:1^D&) :: type_send_r_sts_1
113 integer,
dimension( 0:3^D&) :: type_recv_r_sts_1
114 integer,
dimension( 0:3^D&) :: type_recv_p_sts_1, type_send_p_sts_1
116 integer,
dimension(-1:1^D&) :: type_send_srl_sts_2, type_recv_srl_sts_2
117 integer,
dimension(-1:1^D&) :: type_send_r_sts_2
118 integer,
dimension( 0:3^D&) :: type_recv_r_sts_2
119 integer,
dimension( 0:3^D&) :: type_recv_p_sts_2, type_send_p_sts_2
127 logical :: types_initialized
128 logical :: evolve_magnetic_field
129 procedure(subr1),
pointer,
nopass :: sts_set_sources
130 procedure(subr2),
pointer,
nopass :: sts_getdt
131 procedure(subr5),
pointer,
nopass :: sts_before_first_cycle, sts_after_last_cycle
132 procedure(subr_e),
pointer,
nopass :: sts_handle_errors
133 type(sts_term),
pointer :: next
141 procedure(subr4),
pointer :: sts_get_ncycles
149 if(.not. sts_initialized)
then
152 sts_dtpar=sts_dtpar/dble(
ndim)
153 sts_initialized = .true.
154 if(sts_method .eq. 1)
then
156 sts_get_ncycles => sts_get_ncycles1
157 else if(sts_method .eq. 2)
then
159 sts_get_ncycles => sts_get_ncycles2
161 call mpistop(
"Unknown sts method")
169 if (sts_initialized)
then
177 subroutine sts_params_read(files)
179 character(len=*),
intent(in) :: files(:)
182 namelist /sts_list/ sts_dtpar,sts_ncycles,sts_method,
sourcetype_sts
184 do n = 1,
size(files)
185 open(
unitpar, file=trim(files(n)), status=
"old")
186 read(
unitpar, sts_list,
end=111)
190 end subroutine sts_params_read
199 subroutine add_sts_method(sts_getdt, sts_set_sources, startVar, nflux, startwbc, nwbc, evolve_B)
203 integer,
intent(in) :: startvar, nflux, startwbc, nwbc
204 logical,
intent(in) :: evolve_b
208 subroutine sts_set_sources(ixI^L,ixO^L,w,x,wres,fix_conserve_at_step,my_dt,igrid,nflux)
211 integer,
intent(in) :: ixi^
l, ixo^
l, igrid, nflux
212 double precision,
intent(in) :: x(ixi^s,1:
ndim)
213 double precision,
intent(inout) :: wres(ixi^s,1:nw), w(ixi^s,1:nw)
214 double precision,
intent(in) :: my_dt
215 logical,
intent(in) :: fix_conserve_at_step
216 end subroutine sts_set_sources
218 function sts_getdt(w,ixG^L,ix^L,dx^D,x)
result(dtnew)
220 integer,
intent(in) :: ixg^
l, ix^
l
221 double precision,
intent(in) ::
dx^
d, x(ixg^s,1:
ndim)
222 double precision,
intent(in) :: w(ixg^s,1:nw)
223 double precision :: dtnew
224 end function sts_getdt
228 type(sts_term),
pointer :: temp
231 temp%sts_getdt => sts_getdt
232 temp%sts_set_sources => sts_set_sources
233 temp%sts_before_first_cycle => null()
234 temp%sts_after_last_cycle => null()
235 temp%sts_handle_errors => null()
236 temp%startVar = startvar
237 temp%endVar= startvar+nflux-1
239 temp%startwbc = startwbc
241 temp%types_initialized = .false.
242 temp%evolve_magnetic_field=evolve_b
254 subroutine sts_before_first_cycle(ixI^L, ixO^L, w, x)
256 integer,
intent(in) :: ixi^
l, ixo^
l
257 double precision,
intent(in) :: x(ixi^s,1:
ndim)
258 double precision,
intent(inout) :: w(ixi^s,1:nw)
259 end subroutine sts_before_first_cycle
261 subroutine sts_after_last_cycle(ixI^L, ixO^L, w, x)
263 integer,
intent(in) :: ixi^
l, ixo^
l
264 double precision,
intent(in) :: x(ixi^s,1:
ndim)
265 double precision,
intent(inout) :: w(ixi^s,1:nw)
266 end subroutine sts_after_last_cycle
279 subroutine sts_error_handling(w, x, ixI^L, ixO^L, step)
282 integer,
intent(in) :: ixi^
l,ixo^
l
283 double precision,
intent(inout) :: w(ixi^s,1:nw)
284 double precision,
intent(in) :: x(ixi^s,1:
ndim)
285 integer,
intent(in) :: step
286 end subroutine sts_error_handling
293 function sts_get_ncycles1(dt,dtnew,dt_modified)
result(is)
294 double precision,
intent(in) :: dtnew
295 double precision,
intent(inout) :: dt
296 logical,
intent(inout) :: dt_modified
299 double precision :: ss
302 ss = dtnew*((2.d0*sts_ncycles+1)**2-9.d0)/16.d0
310 if(ss .le. 1.d0)
then
313 is=ceiling((dsqrt(9.d0+16.d0*ss)-1.d0)*0.5d0)
318 end function sts_get_ncycles1
321 function sts_get_ncycles2(dt,dtnew,dt_modified)
result(is)
322 double precision,
intent(in) :: dtnew
323 double precision,
intent(inout) :: dt
324 logical,
intent(inout) :: dt_modified
327 double precision :: ss,rr
332 ncycles = sts_ncycles
334 ss=sum_chev(nu_sts,ncycles,rr)
343 end function sts_get_ncycles2
349 double precision,
intent(inout) :: my_dt
350 double precision :: my_dt1
351 logical :: dt_modified, dt_modified1, dt_modified2
353 double precision :: dtnew,dtmin_mype
354 double precision ::
dx^
d, ss
355 integer:: iigrid, igrid, ncycles
356 type(sts_term),
pointer :: temp,oldtemp
359 dt_modified = .false.
360 do while(
associated(temp))
361 dt_modified2 = .false.
364 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
371 dtmin_mype=min(dtmin_mype, sts_dtpar * temp%sts_getdt(ps(igrid)%w,ixg^
ll,
ixm^
ll,
dx^
d,ps(igrid)%x))
374 call mpi_allreduce(dtmin_mype,dtnew,1,mpi_double_precision,mpi_min,
icomm,
ierrmpi)
375 temp%s = sts_get_ncycles(my_dt,dtnew,dt_modified2)
383 if(dt_modified2)
then
388 dt_modified1 = .false.
389 do while(.not.
associated(oldtemp,temp))
390 oldtemp%s = sts_get_ncycles(my_dt1,oldtemp%dt_expl,dt_modified1)
392 if(dt_modified1)
call mpistop(
"sts dt modified twice")
393 oldtemp=>oldtemp%next
402 pure FUNCTION chev(j,nu,N)
405 double precision,
INTENT(IN) :: nu
406 INTEGER,
INTENT(IN) :: j, n
407 double precision :: chev
409 chev = 1d0 / ((-1d0 + nu)*cos(((2d0*j - 1d0) / n)* (
dpi/2d0)) + 1d0 + nu)
413 FUNCTION sum_chev(nu,N,limMax)
414 double precision,
intent(in) :: nu,limmax
415 integer,
intent(inout) :: n
416 double precision :: sum_chev, tmp
422 do while (j < n .and. sum_chev < limmax)
423 sum_chev = sum_chev + chev(j,nu,n)
427 END FUNCTION sum_chev
442 subroutine sts_add_source2(my_dt)
449 double precision,
intent(in) :: my_dt
450 double precision,
allocatable :: bj(:)
451 double precision :: sumbj,dtj
453 integer:: iigrid, igrid, j, ixc^
l
454 logical :: stagger_flag=.false., prolong_flag=.false., coarsen_flag=.false.
455 type(sts_term),
pointer :: temp
463 do while(
associated(temp))
465 if(.not.temp%evolve_magnetic_field)
then
476 if(
associated(temp%sts_before_first_cycle))
then
481 do iigrid=1,igridstail; igrid=igrids(iigrid);
484 call temp%sts_before_first_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
488 allocate(bj(1:temp%s))
490 bj(j) = chev(j,nu_sts,sts_ncycles)
500 if(.not. temp%types_initialized)
then
502 if(temp%nflux>temp%nwbc)
then
518 temp%types_initialized = .true.
523 if(j .eq. temp%s .and. (sumbj + bj(j)) * temp%dt_expl > my_dt)
then
524 dtj = my_dt - sumbj * temp%dt_expl
526 dtj = bj(j)* temp%dt_expl
528 sumbj = sumbj + bj(j)
531 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
534 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)
535 if(temp%nflux>
ndir)
then
536 ps(igrid)%w(
ixm^t,temp%startVar)=ps(igrid)%w(
ixm^t,temp%startVar)+dtj*ps1(igrid)%w(
ixm^t,temp%startVar)
538 ps(igrid)%ws(ixc^s,1:nws)=ps(igrid)%ws(ixc^s,1:nws)+dtj*ps1(igrid)%w(ixc^s,iw_mag(1:nws))
544 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
547 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)
548 ps(igrid)%w(
ixm^t,temp%startVar:temp%endVar)=ps(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+&
549 dtj*ps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
554 if(fix_conserve_at_step)
then
562 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
568 if(
associated(temp%sts_handle_errors))
then
570 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
573 call temp%sts_handle_errors(ps(igrid)%w,ps(igrid)%x,ixg^
ll,
ixm^
ll,j)
578 if(temp%nflux>temp%nwbc.and.temp%s==j)
then
592 if(
associated(temp%sts_after_last_cycle))
then
593 do iigrid=1,igridstail; igrid=igrids(iigrid);
596 call temp%sts_after_last_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
603 if(.not.temp%evolve_magnetic_field)
then
626 do iigrid=1,igridstail; igrid=igrids(iigrid);
632 end subroutine sts_add_source2
636 subroutine sts_add_source1(my_dt)
644 double precision,
intent(in) :: my_dt
645 double precision :: dtj
646 double precision :: omega1,cmu,cmut,cnu,cnut,one_mu_nu
647 double precision,
allocatable :: bj(:)
648 integer:: iigrid, igrid, j, ixc^
l, ixgext^
l
649 logical :: evenstep, stagger_flag=.false., prolong_flag=.false., coarsen_flag=.false., total_energy_flag=.true.
650 type(sts_term),
pointer :: temp
651 type(state),
dimension(:),
pointer :: tmpps1, tmpps2
659 do while(
associated(temp))
661 if(.not.temp%evolve_magnetic_field)
then
672 if(
associated(temp%sts_before_first_cycle))
then
680 do iigrid=1,igridstail; igrid=igrids(iigrid);
683 call temp%sts_before_first_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
684 if(.not.
allocated(ps2(igrid)%w))
allocate(ps2(igrid)%w(ixg^t,1:nw))
685 if(.not.
allocated(ps3(igrid)%w))
allocate(ps3(igrid)%w(ixg^t,1:nw))
686 if(.not.
allocated(ps4(igrid)%w))
allocate(ps4(igrid)%w(ixg^t,1:nw))
687 ps1(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
688 ps2(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
693 ixgext^
l=ixg^
ll^ladd1;
695 do iigrid=1,igridstail; igrid=igrids(iigrid);
696 if(.not.
allocated(ps2(igrid)%w))
then
699 if(.not.
allocated(ps3(igrid)%w))
allocate(ps3(igrid)%w(ixg^t,1:nw))
700 if(.not.
allocated(ps4(igrid)%w))
allocate(ps4(igrid)%w(ixg^t,1:nw))
701 ps1(igrid)%w=ps(igrid)%w
702 ps2(igrid)%w=ps(igrid)%w
703 ps1(igrid)%ws=ps(igrid)%ws
704 ps2(igrid)%ws=ps(igrid)%ws
709 do iigrid=1,igridstail; igrid=igrids(iigrid);
710 if(.not.
allocated(ps2(igrid)%w))
allocate(ps2(igrid)%w(ixg^t,1:nw))
711 if(.not.
allocated(ps3(igrid)%w))
allocate(ps3(igrid)%w(ixg^t,1:nw))
712 if(.not.
allocated(ps4(igrid)%w))
allocate(ps4(igrid)%w(ixg^t,1:nw))
713 ps1(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
714 ps2(igrid)%w(ixg^t,1:nw)=ps(igrid)%w(ixg^t,1:nw)
720 allocate(bj(0:temp%s))
724 omega1=4.d0/dble(temp%s**2+temp%s-2)
738 if(.not. temp%types_initialized)
then
740 if(temp%nflux>temp%nwbc)
then
756 temp%types_initialized = .true.
761 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
765 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)
767 ps3(igrid)%w(ixc^s,temp%startVar:temp%endVar) = my_dt * ps4(igrid)%w(ixc^s,temp%startVar:temp%endVar)
768 if(temp%nflux>
ndir)
then
769 ps1(igrid)%w(
ixm^t,temp%startVar) = ps1(igrid)%w(
ixm^t,temp%startVar) + cmut * ps3(igrid)%w(
ixm^t,temp%startVar)
771 ps1(igrid)%ws(ixc^s,1:nws) = ps1(igrid)%ws(ixc^s,1:nws) + cmut * ps3(igrid)%w(ixc^s,iw_mag(1:nws))
777 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
780 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)
782 ps3(igrid)%w(
ixm^t,temp%startVar:temp%endVar) = my_dt * ps4(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
783 ps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar) = ps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar) + &
784 cmut * ps3(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
788 if(fix_conserve_at_step)
then
796 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
803 if(
associated(temp%sts_handle_errors))
then
805 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
808 call temp%sts_handle_errors(ps1(igrid)%w,ps1(igrid)%x,ixg^
ll,
ixm^
ll,1)
812 if(temp%nflux>temp%nwbc.and.temp%s==1)
then
831 bj(j)=dble(j**2+j-2)/dble(2*j*(j+1))
832 cmu=dble(2*j-1)/dble(j)*bj(j)/bj(j-1)
834 cnu=dble(1-j)/dble(j)*bj(j)/bj(j-2)
835 cnut=(bj(j-1)-1.d0)*cmut
836 one_mu_nu=1.d0-cmu-cnu
848 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
854 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)
855 if(temp%nflux>
ndir)
then
856 tmpps2(igrid)%w(
ixm^t,temp%startVar)=cmu*tmpps1(igrid)%w(
ixm^t,temp%startVar)+&
857 cnu*tmpps2(igrid)%w(
ixm^t,temp%startVar)+one_mu_nu*ps(igrid)%w(
ixm^t,temp%startVar)+&
858 dtj*ps4(igrid)%w(
ixm^t,temp%startVar)+cnut*ps3(igrid)%w(
ixm^t,temp%startVar)
860 tmpps2(igrid)%ws(ixc^s,1:nws)=cmu*tmpps1(igrid)%ws(ixc^s,1:nws)+&
861 cnu*tmpps2(igrid)%ws(ixc^s,1:nws)+one_mu_nu*ps(igrid)%ws(ixc^s,1:nws)+&
862 dtj*ps4(igrid)%w(ixc^s,iw_mag(1:nws))+cnut*ps3(igrid)%w(ixc^s,iw_mag(1:nws))
868 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
874 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)
875 tmpps2(igrid)%w(
ixm^t,temp%startVar:temp%endVar)=cmu*tmpps1(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+&
876 cnu*tmpps2(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+one_mu_nu*ps(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+&
877 dtj*ps4(igrid)%w(
ixm^t,temp%startVar:temp%endVar)+cnut*ps3(igrid)%w(
ixm^t,temp%startVar:temp%endVar)
881 if(fix_conserve_at_step)
then
889 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
895 if(
associated(temp%sts_handle_errors))
then
897 do iigrid=1,igridstail_active; igrid=igrids_active(iigrid);
900 call temp%sts_handle_errors(tmpps2(igrid)%w,ps(igrid)%x,ixg^
ll,
ixm^
ll,j)
904 if(temp%nflux>temp%nwbc.and.temp%s==j)
then
916 evenstep=.not.evenstep
919 if(
associated(temp%sts_after_last_cycle))
then
921 do iigrid=1,igridstail; igrid=igrids(iigrid);
924 ps(igrid)%w(ixg^t,temp%startVar:temp%endVar)=tmpps2(igrid)%w(ixg^t,temp%startVar:temp%endVar)
925 call temp%sts_after_last_cycle(ixg^
ll,ixg^
ll,ps(igrid)%w,ps(igrid)%x)
934 do iigrid=1,igridstail; igrid=igrids(iigrid);
935 ps(igrid)%w(ixg^t,temp%startVar:temp%endVar)=tmpps2(igrid)%w(ixg^t,temp%startVar:temp%endVar)
936 ps(igrid)%ws=tmpps2(igrid)%ws
941 do iigrid=1,igridstail; igrid=igrids(iigrid);
942 ps(igrid)%w(ixg^t,temp%startVar:temp%endVar)=tmpps2(igrid)%w(ixg^t,temp%startVar:temp%endVar)
950 if(.not.temp%evolve_magnetic_field)
then
973 do iigrid=1,igridstail; igrid=igrids(iigrid);
979 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
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
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