mod_afld Module Reference

Module for including anisotropic flux limited diffusion (AFLD)-approximation in Radiation-hydrodynamics simulations Based on Turner and stone 2001. See [1]Moens, N., Sundqvist, J. O., El Mellah, I., Poniatowski, L., Teunissen, J., and Keppens, R., Radiation-hydrodynamics with MPI-AMRVAC . Flux-limited diffusion Astronomy and Astrophysics, vol. 657, 2022. doi:10.1051/0004-6361/202141023. For more information. with RK updates on 16/03/26. More...

Functions/Subroutines
subroutine	afld_params_read (files)
	public methods these are called in mod_hd_phys or mod_mhd_phys
subroutine, public	afld_init (he_abundance, afld_gamma)
	Initialising FLD-module: Read opacities Initialise Multigrid adimensionalise kappa Add extra variables to w-array, flux, kappa, eddington Tensor Lambda and R ...
subroutine, public	add_afld_rad_force (qdt, ixil, ixol, wct, wctprim, w, x, energy, qsourcesplit, active)
	w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the radiation force
subroutine, public	afld_radforce_get_dt (w, ixil, ixol, dtnew, dxd, x)
	get dt limit for radiation force: NOTE: only uniform cartesian here!
subroutine, public	get_afld_energy_interact (qdt, ixil, ixol, wct, w, x, energy, qsourcesplit, active)
	w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the energy exchange between gas and radiation
subroutine, public	afld_get_opacity (w, x, ixil, ixol, fld_kappa)
	Sets the opacity in the w-array by calling mod_opal_opacity.
subroutine, public	afld_get_fluxlimiter (w, x, ixil, ixol, fld_lambda, fld_r, nth)
	Calculate fld flux limiter This subroutine calculates flux limiter lambda using the prescription stored in fld_fluxlimiter. It also calculates the ratio of radiation scaleheight and mean free path.
subroutine, public	afld_get_radflux (w, x, ixil, ixol, rad_flux, nth)
	Calculate Radiation Flux Calculate the Flux using the fld closure relation F = -c*lambda/(kappa*rho) *grad E.
subroutine	afld_get_eddington (w, x, ixil, ixol, eddington_tensor, nth)
	Calculate Eddington-tensor Stores Eddington-tensor in w-array.
subroutine, public	afld_get_radpress (w, x, ixil, ixol, rad_pressure, nth)
	Calculate Radiation Pressure Returns Radiation Pressure as tensor.
subroutine	diffuse_e_rad_mg (dtfactor, qdt, qtc, psa, psb)
	Calling all subroutines to perform the multigrid method Communicates rad_e and diff_coeff to multigrid library.
subroutine	evaluate_e_rad_mg (qtc, psa)
	inplace update of psa==>F_im(psa)
subroutine	put_diffterm_onegrid (ixil, ixol, w)
	inplace update of psa==>F_im(psa)
subroutine	afld_get_diffcoef_central (w, wct, wctprim, x, ixil, ixol, primitives_filled)
	Calculates cell-centered diffusion coefficient to be used in multigrid.
subroutine	update_diffcoeff (psa)
subroutine	energy_interaction (w, wct, x, ixil, ixol)
	This subroutine calculates the radiation heating, radiation cooling and photon tiring using an implicit scheme. These sourceterms are applied using the root-finding of a 4th order polynomial This routine loops over every cell in the domain and computes the coefficients of the polynomials in every cell.
subroutine	bisection_method (e_gas, e_rad, c0, c1)
	Find the root of the 4th degree polynomial using the bisection method.
subroutine	newton_method (e_gas, e_rad, c0, c1)
	Find the root of the 4th degree polynomial using the Newton-Ralphson method.
subroutine	halley_method (e_gas, e_rad, c0, c1)
	Find the root of the 4th degree polynomial using the Halley method.
double precision function	polynomial_bisection (e_gas, c0, c1)
	Evaluate polynomial at argument e_gas.
double precision function	dpolynomial_bisection (e_gas, c0, c1)
	Evaluate first derivative of polynomial at argument e_gas.
double precision function	ddpolynomial_bisection (e_gas, c0, c1)
	Evaluate second derivative of polynomial at argument e_gas.

Variables
logical	fld_eint_split = .false.
	source split for energy interact and radforce:
logical	fld_radforce_split = .false.
double precision, public	fld_kappa0 = 0.34d0
	Opacity per unit of unit_density.
double precision, public	afld_mu = 0.6d0
	mean particle mass
double precision, public	fld_bisect_tol = 1.d-4
	Tolerance for bisection method for Energy sourceterms This is a percentage of the minimum of gas- and radiation energy.
double precision, public	fld_diff_tol = 1.d-4
	Tolerance for adi method for radiative Energy diffusion.
character(len=8), dimension(:), allocatable	fld_opacity_law
	opacity switches
character(len=50)	fld_opal_table = 'Y09800'
character(len=16)	fld_fluxlimiter = 'Pomraning'
	flux limiter switch
integer, dimension(:), allocatable	i_diff_mg
	diffusion coefficient for multigrid method
character(len=8)	fld_interaction_method = 'Halley'
	Which method to find the root for the energy interaction polynomial.
logical	diffcrash_resume = .true.
	Resume run when multigrid returns error.
double precision	dt_diff = 0.d0
	running timestep for diffusion solver, initialised as zero

Detailed Description

Module for including anisotropic flux limited diffusion (AFLD)-approximation in Radiation-hydrodynamics simulations Based on Turner and stone 2001. See [1]Moens, N., Sundqvist, J. O., El Mellah, I., Poniatowski, L., Teunissen, J., and Keppens, R., Radiation-hydrodynamics with MPI-AMRVAC . Flux-limited diffusion Astronomy and Astrophysics, vol. 657, 2022. doi:10.1051/0004-6361/202141023. For more information. with RK updates on 16/03/26.

Function/Subroutine Documentation

add_afld_rad_force()

subroutine, public mod_afld::add_afld_rad_force	(	double precision, intent(in)	qdt,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixi^s,1:nw), intent(in)	wct,
		double precision, dimension(ixi^s,1:nw), intent(in)	wctprim,
		double precision, dimension(ixi^s,1:nw), intent(inout)	w,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x,
		logical, intent(in)	energy,
		logical, intent(in)	qsourcesplit,
		logical, intent(inout)	active
	)

w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the radiation force

Calculate and add sourceterms

Radiation force = kappa*rho/c *Flux = lambda gradE

Momentum equation source term

Energy equation source term (kinetic energy)

Photon tiring calculate tensor div_v !$OMP PARALLEL DO

!$OMP END PARALLEL DO

VARIABLE NAMES DIV ARE ACTUALLY GRADIENTS

eq 34 Turner and stone (Only 2D)

Definition at line 116 of file mod_afld.t.

Here is the call graph for this function:

afld_get_diffcoef_central()

subroutine mod_afld::afld_get_diffcoef_central	(	double precision, dimension(ixi^s, 1:nw), intent(inout)	w,
		double precision, dimension(ixi^s, 1:nw), intent(in)	wct,
		double precision, dimension(ixi^s, 1:nw), intent(in)	wctprim,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		logical, intent(in)	primitives_filled
	)

Calculates cell-centered diffusion coefficient to be used in multigrid.

calculate diffusion coefficient

Definition at line 607 of file mod_afld.t.

Here is the call graph for this function:

afld_get_eddington()

subroutine mod_afld::afld_get_eddington	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim,1:ndim), intent(out)	eddington_tensor,
		integer, intent(in)	nth
	)

Calculate Eddington-tensor Stores Eddington-tensor in w-array.

Calculate R everywhere |grad E|/(rho kappa E)

Calculate radiation pressure P = (lambda + lambda^2 R^2)*E

Definition at line 389 of file mod_afld.t.

Here is the call graph for this function:

afld_get_fluxlimiter()

subroutine, public mod_afld::afld_get_fluxlimiter	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim), intent(out)	fld_lambda,
		double precision, dimension(ixo^s,1:ndim), intent(out)	fld_r,
		integer, intent(in)	nth
	)

Calculate fld flux limiter This subroutine calculates flux limiter lambda using the prescription stored in fld_fluxlimiter. It also calculates the ratio of radiation scaleheight and mean free path.

Calculate R everywhere |grad E|/(rho kappa E)

Calculate the flux limiter, lambda

Calculate R everywhere |grad E|/(rho kappa E)

Calculate the flux limiter, lambda Levermore and Pomraning: lambda = (2 + R)/(6 + 3R + R^2)

Calculate R everywhere |grad E|/(rho kappa E)

Calculate the flux limiter, lambda Minerbo:

Definition at line 282 of file mod_afld.t.

Here is the call graph for this function:

afld_get_opacity()

subroutine, public mod_afld::afld_get_opacity	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim), intent(out)	fld_kappa
	)

Sets the opacity in the w-array by calling mod_opal_opacity.

Definition at line 242 of file mod_afld.t.

Here is the call graph for this function:

afld_get_radflux()

subroutine, public mod_afld::afld_get_radflux	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s, 1:ndim), intent(out)	rad_flux,
		integer, intent(in)	nth
	)

Calculate Radiation Flux Calculate the Flux using the fld closure relation F = -c*lambda/(kappa*rho) *grad E.

Definition at line 364 of file mod_afld.t.

Here is the call graph for this function:

afld_get_radpress()

subroutine, public mod_afld::afld_get_radpress	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim,1:ndim), intent(out)	rad_pressure,
		integer, intent(in)	nth
	)

Calculate Radiation Pressure Returns Radiation Pressure as tensor.

Definition at line 442 of file mod_afld.t.

Here is the call graph for this function:

afld_init()

subroutine, public mod_afld::afld_init	(	double precision, intent(in)	he_abundance,
		double precision, intent(in)	afld_gamma
	)

Initialising FLD-module: Read opacities Initialise Multigrid adimensionalise kappa Add extra variables to w-array, flux, kappa, eddington Tensor Lambda and R ...

allocate boundary conditions

read par files

Need mean molecular weight

Definition at line 79 of file mod_afld.t.

Here is the call graph for this function:

afld_params_read()

subroutine mod_afld::afld_params_read

(

character(len=*), dimension(:), intent(in)

files

)

public methods these are called in mod_hd_phys or mod_mhd_phys

Reading in fld-list parameters from .par file

Definition at line 54 of file mod_afld.t.

afld_radforce_get_dt()

subroutine, public mod_afld::afld_radforce_get_dt	(	double precision, dimension(ixi^s,1:nw), intent(in)	w,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, intent(inout)	dtnew,
		double precision, intent(in)	dx,
		double precision, intent(in)	d,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x
	)

get dt limit for radiation force: NOTE: only uniform cartesian here!

Definition at line 194 of file mod_afld.t.

Here is the call graph for this function:

bisection_method()

subroutine mod_afld::bisection_method	(	double precision, intent(inout)	e_gas,
		double precision, intent(in)	e_rad,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Find the root of the 4th degree polynomial using the bisection method.

Definition at line 734 of file mod_afld.t.

Here is the call graph for this function:

ddpolynomial_bisection()

double precision function mod_afld::ddpolynomial_bisection	(	double precision, intent(in)	e_gas,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Evaluate second derivative of polynomial at argument e_gas.

Definition at line 876 of file mod_afld.t.

diffuse_e_rad_mg()

subroutine mod_afld::diffuse_e_rad_mg	(	double precision, intent(in)	dtfactor,
		double precision, intent(in)	qdt,
		double precision, intent(in)	qtc,
		type(state), dimension(max_blocks), target	psa,
		type(state), dimension(max_blocks), target	psb
	)

Calling all subroutines to perform the multigrid method Communicates rad_e and diff_coeff to multigrid library.

Set diffusion timestep, add previous timestep if mg did not converge:

replace call set_rhs(mg, -1/dt, 0.0_dp)

Reset dt_diff when diffusion worked out

Definition at line 461 of file mod_afld.t.

Here is the call graph for this function:

dpolynomial_bisection()

double precision function mod_afld::dpolynomial_bisection	(	double precision, intent(in)	e_gas,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Evaluate first derivative of polynomial at argument e_gas.

Definition at line 866 of file mod_afld.t.

energy_interaction()

subroutine mod_afld::energy_interaction	(	double precision, dimension(ixi^s,1:nw), intent(inout)	w,
		double precision, dimension(ixi^s,1:nw), intent(in)	wct,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l
	)

This subroutine calculates the radiation heating, radiation cooling and photon tiring using an implicit scheme. These sourceterms are applied using the root-finding of a 4th order polynomial This routine loops over every cell in the domain and computes the coefficients of the polynomials in every cell.

e_gas is the INTERNAL ENERGY without KINETIC ENERGY

Calculate coefficients for polynomial

Loop over every cell for rootfinding method

advance E_rad

Update gas-energy in w, internal + kinetic

Beginning of module substracted wCT Ekin So now add wCT Ekin

Update rad-energy in w

Definition at line 655 of file mod_afld.t.

Here is the call graph for this function:

evaluate_e_rad_mg()

subroutine mod_afld::evaluate_e_rad_mg	(	double precision, intent(in)	qtc,
		type(state), dimension(max_blocks), target	psa
	)

inplace update of psa==>F_im(psa)

Definition at line 565 of file mod_afld.t.

Here is the call graph for this function:

get_afld_energy_interact()

subroutine, public mod_afld::get_afld_energy_interact	(	double precision, intent(in)	qdt,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixi^s,1:nw), intent(in)	wct,
		double precision, dimension(ixi^s,1:nw), intent(inout)	w,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x,
		logical, intent(in)	energy,
		logical, intent(in)	qsourcesplit,
		logical, intent(inout)	active
	)

w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the energy exchange between gas and radiation

Calculate and add sourceterms

Add energy sourceterms

Definition at line 220 of file mod_afld.t.

Here is the call graph for this function:

halley_method()

subroutine mod_afld::halley_method	(	double precision, intent(inout)	e_gas,
		double precision, intent(in)	e_rad,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Find the root of the 4th degree polynomial using the Halley method.

Compare error with dx = dx/dy dy

Definition at line 825 of file mod_afld.t.

Here is the call graph for this function:

newton_method()

subroutine mod_afld::newton_method	(	double precision, intent(inout)	e_gas,
		double precision, intent(in)	e_rad,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Find the root of the 4th degree polynomial using the Newton-Ralphson method.

Compare error with dx = dx/dy dy

Definition at line 796 of file mod_afld.t.

Here is the call graph for this function:

polynomial_bisection()

double precision function mod_afld::polynomial_bisection	(	double precision, intent(in)	e_gas,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Evaluate polynomial at argument e_gas.

Definition at line 856 of file mod_afld.t.

put_diffterm_onegrid()

subroutine mod_afld::put_diffterm_onegrid	(	integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixi^s, 1:nw), intent(inout)	w
	)

inplace update of psa==>F_im(psa)

Definition at line 586 of file mod_afld.t.

update_diffcoeff()

subroutine mod_afld::update_diffcoeff

(

type(state), dimension(max_blocks), target

psa

)

Definition at line 635 of file mod_afld.t.

Here is the call graph for this function:

Variable Documentation

afld_mu

double precision, public mod_afld::afld_mu = 0.6d0

mean particle mass

Definition at line 19 of file mod_afld.t.

diffcrash_resume

logical mod_afld::diffcrash_resume = .true.

Resume run when multigrid returns error.

Definition at line 35 of file mod_afld.t.

dt_diff

double precision mod_afld::dt_diff = 0.d0

running timestep for diffusion solver, initialised as zero

Definition at line 39 of file mod_afld.t.

fld_bisect_tol

double precision, public mod_afld::fld_bisect_tol = 1.d-4

Tolerance for bisection method for Energy sourceterms This is a percentage of the minimum of gas- and radiation energy.

Definition at line 22 of file mod_afld.t.

fld_diff_tol

double precision, public mod_afld::fld_diff_tol = 1.d-4

Tolerance for adi method for radiative Energy diffusion.

Definition at line 24 of file mod_afld.t.

fld_eint_split

logical mod_afld::fld_eint_split = .false.

source split for energy interact and radforce:

Definition at line 14 of file mod_afld.t.

fld_fluxlimiter

character(len=16) mod_afld::fld_fluxlimiter = 'Pomraning'

flux limiter switch

Definition at line 29 of file mod_afld.t.

fld_interaction_method

character(len=8) mod_afld::fld_interaction_method = 'Halley'

Which method to find the root for the energy interaction polynomial.

Definition at line 33 of file mod_afld.t.

fld_kappa0

double precision, public mod_afld::fld_kappa0 = 0.34d0

Opacity per unit of unit_density.

Definition at line 17 of file mod_afld.t.

fld_opacity_law

character(len=8), dimension(:), allocatable mod_afld::fld_opacity_law

opacity switches

Definition at line 26 of file mod_afld.t.

fld_opal_table

character(len=50) mod_afld::fld_opal_table = 'Y09800'

Definition at line 27 of file mod_afld.t.

fld_radforce_split

logical mod_afld::fld_radforce_split = .false.

Definition at line 15 of file mod_afld.t.

i_diff_mg

integer, dimension(:), allocatable mod_afld::i_diff_mg

diffusion coefficient for multigrid method

Definition at line 31 of file mod_afld.t.