module mms_braginskii_m !! Implementation of Method of Manufactured Solutions for BRAGINSKII model in circular geometry use MPI use comm_handler_m, only : comm_handler_t use precision_grillix_m, only : GP, MPI_GP, GP_EPS use error_handling_grillix_m, only: handle_error use status_codes_grillix_m, only : GRILLIX_ERR_NAMELIST, GRILLIX_ERR_OTHER use variable_m, only : variable_t use parallel_map_m, only : parallel_map_t use commons_misc_m, only : error_ananum use constants_m, only : Pi, TWO_PI use screen_io_m, only : get_stdout use parallelisation_setup_m, only : is_rank_info_writer use circular_equilibrium_m, only : circular_t use slab_equilibrium_m, only : slab_t use dommaschk_equilibrium_m, only : dommaschk_t use equilibrium_m, only : equilibrium_t use mesh_cart_m, only : mesh_cart_t use mms_braginskii_circular_m use mms_braginskii_slab_m use mms_braginskii_dommaschk_m ! Parameters use params_brag_model_m, only : & beta, mass_ratio_ei, boussinesq_on implicit none public :: mms_diagnostics_braginskii public :: mms_sol_braginskii_ne public :: mms_sol_braginskii_pot public :: mms_sol_braginskii_vort public :: mms_sol_braginskii_upar public :: mms_sol_braginskii_jpar public :: mms_sol_braginskii_apar public :: mms_sol_braginskii_psipar public :: mms_sol_braginskii_te public :: mms_sol_braginskii_ti public :: mms_sol_braginskii_qe public :: mms_sol_braginskii_qi public :: mms_source_braginskii_continuity public :: mms_source_braginskii_vorticity public :: mms_source_braginskii_parmomentum public :: mms_source_braginskii_ohm public :: mms_source_braginskii_etemp public :: mms_source_braginskii_itemp public :: mms_source_braginskii_landau_e public :: mms_source_braginskii_landau_i contains subroutine mms_diagnostics_braginskii(comm_handler, equi, mesh_cano, mesh_stag, map, tau, & ne, pot, vort, upar, jpar, apar, te, ti) !! Prints information on numerical errors of numerical solutions type(comm_handler_t), intent(in) :: comm_handler !! Communicators class(equilibrium_t), intent(inout) :: equi !! Equilibrium type(mesh_cart_t), intent(in) :: mesh_cano !! Mesh (canonical) type(mesh_cart_t), intent(in) :: mesh_stag !! Mesh (staggered) type(parallel_map_t), intent(in) :: map !! Mesh within poloidal plane real(GP), intent(in) :: tau !! Time type(variable_t), intent(in) :: ne !! Numerical solution for electron density type(variable_t), intent(in) :: pot !! Numerical solution for electrostatic potential type(variable_t), intent(in) :: vort !! Numerical solution for generalised vorticity type(variable_t), intent(in) :: upar !! Numerical solution for parallel ion velocity type(variable_t), intent(in) :: jpar !! Numerical solution for parallel current type(variable_t), intent(in) :: apar !! Numerical solution for parallel electromagnetic potential type(variable_t), intent(in) :: te !! Numerical solution for electron temperature type(variable_t), intent(in) :: ti !! Numerical solution for ion temperature integer :: l real(GP) :: x, y, phi_cano, phi_stag real(GP), dimension(mesh_cano%get_n_points()) :: ne_mms, pot_mms, vort_mms, te_mms, ti_mms real(GP), dimension(mesh_stag%get_n_points()) :: upar_mms, jpar_mms, apar_mms real(GP) :: ne_nrm2, ne_nrmsup, ne_err2, ne_errsup real(GP) :: pot_nrm2, pot_nrmsup, pot_err2, pot_errsup real(GP) :: vort_nrm2, vort_nrmsup, vort_err2, vort_errsup real(GP) :: upar_nrm2, upar_nrmsup, upar_err2, upar_errsup real(GP) :: jpar_nrm2, jpar_nrmsup, jpar_err2, jpar_errsup real(GP) :: apar_nrm2, apar_nrmsup, apar_err2, apar_errsup real(GP) :: te_nrm2, te_nrmsup, te_err2, te_errsup real(GP) :: ti_nrm2, ti_nrmsup, ti_err2, ti_errsup character(len=22) :: fname logical :: fexist integer :: ioerr character(len=256) :: io_errmsg phi_cano = mesh_cano%get_phi() phi_stag = mesh_stag%get_phi() ! Get MMS solutions ----------------------------------------------------------------------- !$omp parallel default(none) private(l, x, y) & !$omp shared(equi, tau, mesh_cano, phi_cano, mesh_stag, phi_stag, & !$omp ne_mms, pot_mms, vort_mms, te_mms, ti_mms, & !$omp upar_mms, jpar_mms, apar_mms) !$omp do do l = 1, mesh_cano%get_n_points() x = mesh_cano%get_x(l) y = mesh_cano%get_y(l) ne_mms(l) = mms_sol_braginskii_ne(equi, x, y, phi_cano, tau) pot_mms(l) = mms_sol_braginskii_pot(equi, x, y, phi_cano, tau) vort_mms(l) = mms_sol_braginskii_vort(equi, x, y, phi_cano, tau) te_mms(l) = mms_sol_braginskii_te(equi, x, y, phi_cano, tau) ti_mms(l) = mms_sol_braginskii_ti(equi, x, y, phi_cano, tau) enddo !$omp end do !$omp do do l = 1, mesh_stag%get_n_points() x = mesh_stag%get_x(l) y = mesh_stag%get_y(l) upar_mms(l) = mms_sol_braginskii_upar(equi, x, y, phi_stag, tau) jpar_mms(l) = mms_sol_braginskii_jpar(equi, x, y, phi_stag, tau) apar_mms(l) = mms_sol_braginskii_apar(equi, x, y, phi_stag, tau) enddo !$omp end do !$omp end parallel ! Compute errors between MMS and numeric solution ----------------------------------------- call error_ananum(comm_handler, mesh_cano%get_n_points(), ne_mms, ne%vals, ne_nrm2, ne_nrmsup, ne_err2, ne_errsup) call error_ananum(comm_handler, mesh_cano%get_n_points(), pot_mms, pot%vals, pot_nrm2, pot_nrmsup, pot_err2, pot_errsup) call error_ananum(comm_handler, mesh_cano%get_n_points(), vort_mms, vort%vals, vort_nrm2, vort_nrmsup, vort_err2, vort_errsup) call error_ananum(comm_handler, mesh_stag%get_n_points(), upar_mms, upar%vals, upar_nrm2, upar_nrmsup, upar_err2, upar_errsup) call error_ananum(comm_handler, mesh_stag%get_n_points(), jpar_mms, jpar%vals, jpar_nrm2, jpar_nrmsup, jpar_err2, jpar_errsup) call error_ananum(comm_handler, mesh_stag%get_n_points(), apar_mms, apar%vals, apar_nrm2, apar_nrmsup, apar_err2, apar_errsup) call error_ananum(comm_handler, mesh_cano%get_n_points(), te_mms, te%vals, te_nrm2, te_nrmsup, te_err2, te_errsup) call error_ananum(comm_handler, mesh_cano%get_n_points(), ti_mms, ti%vals, ti_nrm2, ti_nrmsup, ti_err2, ti_errsup) ! Write out results on display ------------------------------------------------------------ if (is_rank_info_writer) then write(get_stdout(),*)'MMS diagnostics ---------------------------------------' write(get_stdout(),202)tau 202 FORMAT(3X,'tau = ',ES14.6E3) write(get_stdout(),203)'Field','norm(L2)','norm(Linf)','err(L2)','errnrm(L2)', 'err(Linf)', 'errnrm(Linf)' 203 FORMAT(A8,6(A19)) write(get_stdout(),204)'ne: ', ne_nrm2, ne_nrmsup, ne_err2, ne_err2/(ne_nrm2+GP_EPS), & ne_errsup, ne_errsup/(ne_nrmsup+GP_EPS) write(get_stdout(),204)'pot: ', pot_nrm2, pot_nrmsup, pot_err2, pot_err2/(pot_nrm2+GP_EPS), & pot_errsup, pot_errsup/(pot_nrmsup+GP_EPS) write(get_stdout(),204)'vort: ', vort_nrm2, vort_nrmsup, vort_err2, vort_err2/(vort_nrm2+GP_EPS), & vort_errsup, vort_errsup/(vort_nrmsup+GP_EPS) write(get_stdout(),204)'upar: ', upar_nrm2, upar_nrmsup, upar_err2, upar_err2/(upar_nrm2+GP_EPS), & upar_errsup, upar_errsup/(upar_nrmsup+GP_EPS) write(get_stdout(),204)'jpar: ', jpar_nrm2, jpar_nrmsup, jpar_err2, jpar_err2/(jpar_nrm2+GP_EPS), & jpar_errsup, jpar_errsup/(jpar_nrmsup+GP_EPS) write(get_stdout(),204)'apar: ', apar_nrm2, apar_nrmsup, apar_err2, apar_err2/(apar_nrm2+GP_EPS), & apar_errsup, apar_errsup/(apar_nrmsup+GP_EPS) write(get_stdout(),204)'te: ', te_nrm2, te_nrmsup, te_err2, te_err2/(te_nrm2+GP_EPS), & te_errsup, te_errsup/(te_nrmsup+GP_EPS) write(get_stdout(),204)'ti: ', ti_nrm2, ti_nrmsup, ti_err2, ti_err2/(ti_nrm2+GP_EPS), & ti_errsup, ti_errsup/(ti_nrmsup+GP_EPS) 204 FORMAT(3X,A8,6(5X,ES14.6E3)) write(get_stdout(),*)'-------------------------------------------------------' endif ! Write out result to file for integration tests ------------------------------------------ if (is_rank_info_writer) then fname = 'mms_braginskii_job.out' inquire(file=fname, exist = fexist) if (fexist) then open(unit=56, file=fname, status='old', position='append', & action='write', iostat=ioerr, iomsg=io_errmsg) else open(unit=56, file=fname, status='new', & action='write', iostat=ioerr, iomsg=io_errmsg) endif if (ioerr /= 0) then call handle_error(io_errmsg, GRILLIX_ERR_NAMELIST, __LINE__, __FILE__) endif write(56,603)tau write(56,603)ne_nrm2 write(56,603)ne_nrmsup write(56,603)ne_err2 write(56,603)ne_errsup write(56,603)pot_nrm2 write(56,603)pot_nrmsup write(56,603)pot_err2 write(56,603)pot_errsup write(56,603)vort_nrm2 write(56,603)vort_nrmsup write(56,603)vort_err2 write(56,603)vort_errsup write(56,603)upar_nrm2 write(56,603)upar_nrmsup write(56,603)upar_err2 write(56,603)upar_errsup write(56,603)jpar_nrm2 write(56,603)jpar_nrmsup write(56,603)jpar_err2 write(56,603)jpar_errsup write(56,603)apar_nrm2 write(56,603)apar_nrmsup write(56,603)apar_err2 write(56,603)apar_errsup write(56,603)te_nrm2 write(56,603)te_nrmsup write(56,603)te_err2 write(56,603)te_errsup write(56,603)ti_nrm2 write(56,603)ti_nrmsup write(56,603)ti_err2 write(56,603)ti_errsup 603 FORMAT(ES20.12E3) close(56) 301 FORMAT(ES20.12E3) endif end subroutine real(GP) function mms_sol_braginskii_ne(equi, x, y, z, t) !! MMS solution for electron density class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_ne = mms_sol_braginskii_circular_ne(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_ne = mms_sol_braginskii_slab_ne(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_ne = mms_sol_braginskii_dommaschk_ne(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_pot(equi, x, y, z, t) !! MMS solution for electrostatic potential class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_pot = mms_sol_braginskii_circular_pot(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_pot = mms_sol_braginskii_slab_pot(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_pot = mms_sol_braginskii_dommaschk_pot(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_vort(equi, x, y, z, t) !! MMS solution for generalised vorticity class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_vort = mms_sol_braginskii_circular_vort(equi, x, y, z, t) type is(slab_t) if (boussinesq_on) then mms_sol_braginskii_vort = & mms_sol_braginskii_slab_vortbsq(equi, x, y, z, t) else mms_sol_braginskii_vort = & mms_sol_braginskii_slab_vort(equi, x, y, z, t) endif type is(dommaschk_t) mms_sol_braginskii_vort = & mms_sol_braginskii_dommaschk_vort(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_upar(equi, x, y, z, t) !! MMS solution for parallel ion velocity class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_upar = mms_sol_braginskii_circular_upar(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_upar = mms_sol_braginskii_slab_upar(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_upar = mms_sol_braginskii_dommaschk_upar(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_jpar(equi, x, y, z, t) !! MMS solution for parallel current class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_jpar = mms_sol_braginskii_circular_jpar(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_jpar = mms_sol_braginskii_slab_jpar(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_jpar = mms_sol_braginskii_dommaschk_jpar(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_apar(equi, x, y, z, t) !! MMS solution for parallel electromagnetic potential class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_apar = mms_sol_braginskii_circular_apar(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_apar = mms_sol_braginskii_slab_apar(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_apar = mms_sol_braginskii_dommaschk_apar(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_psipar(equi, x, y, z, t) !! MMS solution for generalised electromagnetic potential (aixiliary) class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP) :: mu mu = mass_ratio_ei mms_sol_braginskii_psipar = beta * mms_sol_braginskii_apar(equi, x, y, z, t) & + mu * mms_sol_braginskii_jpar(equi, x, y, z, t) / & mms_sol_braginskii_ne(equi, x, y, z, t) end function real(GP) function mms_sol_braginskii_te(equi, x, y, z, t) !! MMS solution for electron temperature class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_te = mms_sol_braginskii_circular_te(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_te = mms_sol_braginskii_slab_te(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_te = mms_sol_braginskii_dommaschk_te(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_ti(equi, x, y, z, t) !! MMS solution for ion temperature class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) mms_sol_braginskii_ti = mms_sol_braginskii_circular_ti(equi, x, y, z, t) type is(slab_t) mms_sol_braginskii_ti = mms_sol_braginskii_slab_ti(equi, x, y, z, t) type is(dommaschk_t) mms_sol_braginskii_ti = mms_sol_braginskii_dommaschk_ti(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_qe(equi, x, y, z, t) !! MMS solution for electron heat flux class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) call handle_error('Circular equilibrium not available', GRILLIX_ERR_OTHER, __LINE__, __FILE__) type is(slab_t) mms_sol_braginskii_qe = mms_sol_braginskii_slab_qe(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_sol_braginskii_qi(equi, x, y, z, t) !! MMS solution for ion heat flux class(equilibrium_t) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time select type(equi) type is(circular_t) call handle_error('Circular equilibrium not available', GRILLIX_ERR_OTHER, __LINE__, __FILE__) type is(slab_t) mms_sol_braginskii_qi = mms_sol_braginskii_slab_qi(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_continuity(equi, x, y, z, t, chi, lognepen, src_ne) !! MMS source for continuity equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: chi !! Value of characteristc function of penalisation real(GP), intent(in) :: lognepen !! Value that shall be penalised to real(GP), intent(in) :: src_ne !! Particle source value select type(equi) type is(circular_t) mms_source_braginskii_continuity = & mms_source_braginskii_circular_continuity(equi, x, y, z, t, chi, lognepen, src_ne) type is(slab_t) mms_source_braginskii_continuity = & mms_source_braginskii_slab_continuity(equi, x, y, z, t, src_ne) type is(dommaschk_t) mms_source_braginskii_continuity = & mms_source_braginskii_dommaschk_continuity(equi, x, y, z, t, src_ne) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_vorticity(equi, x, y, z, t, chi, vortpen, src_vort) !! MMS source for continuity equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: chi !! Value of characteristc function of penalisation real(GP), intent(in) :: vortpen !! Value that shall be penalised to real(GP), intent(in) :: src_vort !! Vorticity source value select type(equi) type is(circular_t) mms_source_braginskii_vorticity = & mms_source_braginskii_circular_vorticity(equi, x, y, z, t, chi, vortpen, src_vort) type is(slab_t) if (boussinesq_on) then mms_source_braginskii_vorticity = & mms_source_braginskii_slab_vorticitybsq(equi, x, y, z, t, src_vort) else mms_source_braginskii_vorticity = & mms_source_braginskii_slab_vorticity(equi, x, y, z, t, src_vort) endif type is(dommaschk_t) mms_source_braginskii_vorticity = & mms_source_braginskii_dommaschk_vorticity(equi, x, y, z, t, src_vort) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_parmomentum(equi, x, y, z, t, chi, uparpen, src_upar) !! MMS source for parallel momentum equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: chi !! Value of characteristc function of penalisation real(GP), intent(in) :: uparpen !! Value that shall be penalised to real(GP), intent(in) :: src_upar !! Parallel momentum source value select type(equi) type is(circular_t) mms_source_braginskii_parmomentum = & mms_source_braginskii_circular_parmomentum(equi, x, y, z, t, chi, uparpen, src_upar) type is(slab_t) mms_source_braginskii_parmomentum = & mms_source_braginskii_slab_parmomentum(equi, x, y, z, t, src_upar) type is(dommaschk_t) mms_source_braginskii_parmomentum = & mms_source_braginskii_dommaschk_parmomentum(equi, x, y, z, t, src_upar) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_ohm(equi, x, y, z, t, chi, psiparpen) !! MMS source for Ohm's law class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: chi !! Value of characteristc function of penalisation real(GP), intent(in) :: psiparpen !! Value that shall be penalised to select type(equi) type is(circular_t) mms_source_braginskii_ohm = & mms_source_braginskii_circular_ohm(equi, x, y, z, t, chi, psiparpen) type is(slab_t) mms_source_braginskii_ohm = & mms_source_braginskii_slab_ohm(equi, x, y, z, t) type is(dommaschk_t) mms_source_braginskii_ohm = & mms_source_braginskii_dommaschk_ohm(equi, x, y, z, t) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_etemp(equi, x, y, z, t, chi, logtepen, src_te) !! MMS source for electron temperature equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: chi !! Value of characteristc function of penalisation real(GP), intent(in) :: logtepen !! Value that shall be penalised to real(GP), intent(in) :: src_te !! Electron temperature source value select type(equi) type is(circular_t) mms_source_braginskii_etemp = & mms_source_braginskii_circular_etemp(equi, x, y, z, t, chi, logtepen, src_te) type is(slab_t) mms_source_braginskii_etemp = & mms_source_braginskii_slab_etemp(equi, x, y, z, t, src_te) type is(dommaschk_t) mms_source_braginskii_etemp = & mms_source_braginskii_dommaschk_etemp(equi, x, y, z, t, src_te) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_itemp(equi, x, y, z, t, chi, logtipen, src_ti) !! MMS source for ion temperature equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: chi !! Value of characteristc function of penalisation real(GP), intent(in) :: logtipen !! Value that shall be penalised to real(GP), intent(in) :: src_ti !! Ion temperature source value select type(equi) type is(circular_t) mms_source_braginskii_itemp = & mms_source_braginskii_circular_itemp(equi, x, y, z, t, chi, logtipen, src_ti) type is(slab_t) mms_source_braginskii_itemp = & mms_source_braginskii_slab_itemp(equi, x, y, z, t, src_ti) type is(dommaschk_t) mms_source_braginskii_itemp = & mms_source_braginskii_dommaschk_itemp(equi, x, y, z, t, src_ti) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_landau_e(equi, x, y, z, t, landaualpha, landaubeta) !! MMS source for elliptic Landau heat flux equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: landaualpha !! Numerical parameter alpha of current Lorentzian real(GP), intent(in) :: landaubeta !! Numerical parameter beta of current Lorentzian select type(equi) type is(circular_t) call handle_error('Circular equilibrium not available', GRILLIX_ERR_OTHER, __LINE__, __FILE__) type is(slab_t) mms_source_braginskii_landau_e = & mms_source_braginskii_slab_landau_e(equi, x, y, z, t, landaualpha, landaubeta) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function real(GP) function mms_source_braginskii_landau_i(equi, x, y, z, t, landaualpha, landaubeta) !! MMS source for elliptic Landau heat flux equation class(equilibrium_t), intent(inout) :: equi !! Equilibrium real(GP), intent(in) :: x !! x-coordinate real(GP), intent(in) :: y !! y-coordinate real(GP), intent(in) :: z !! Axial coordinate z (represented by toroidal coordinate) real(GP), intent(in) :: t !! Time real(GP), intent(in) :: landaualpha !! Numerical parameter alpha of current Lorentzian real(GP), intent(in) :: landaubeta !! Numerical parameter beta of current Lorentzian select type(equi) type is(circular_t) call handle_error('Circular equilibrium not available', GRILLIX_ERR_OTHER, __LINE__, __FILE__) type is(slab_t) mms_source_braginskii_landau_i = & mms_source_braginskii_slab_landau_i(equi, x, y, z, t, landaualpha, landaubeta) class default call handle_error('Equilibrium not valid', GRILLIX_ERR_OTHER, __LINE__, __FILE__) end select end function end module