read_equilibrium Subroutine

public subroutine read_equilibrium()

Reads in the interpolation grid, plasma parameters, and fields and stores the quantities in inter_grid and equil

Arguments

None

Calls

proc~~read_equilibrium~~CallsGraph proc~read_equilibrium read_equilibrium h5ltread_dataset_int_f h5ltread_dataset_int_f proc~read_equilibrium->h5ltread_dataset_int_f h5close_f h5close_f proc~read_equilibrium->h5close_f h5gopen_f h5gopen_f proc~read_equilibrium->h5gopen_f h5gclose_f h5gclose_f proc~read_equilibrium->h5gclose_f h5open_f h5open_f proc~read_equilibrium->h5open_f h5ltread_dataset_double_f h5ltread_dataset_double_f proc~read_equilibrium->h5ltread_dataset_double_f interface~deriv deriv proc~read_equilibrium->interface~deriv proc~h5ltread_dataset_int_scalar_f h5ltread_dataset_int_scalar_f proc~read_equilibrium->proc~h5ltread_dataset_int_scalar_f h5fclose_f h5fclose_f proc~read_equilibrium->h5fclose_f h5fopen_f h5fopen_f proc~read_equilibrium->h5fopen_f proc~h5ltread_dataset_int_scalar_f->h5ltread_dataset_int_f

Called by

proc~~read_equilibrium~~CalledByGraph proc~read_equilibrium read_equilibrium program~fidasim fidasim program~fidasim->proc~read_equilibrium

Contents

Source Code


Source Code

subroutine read_equilibrium
    !+ Reads in the interpolation grid, plasma parameters, and fields
    !+ and stores the quantities in [[libfida:inter_grid]] and [[libfida:equil]]
    integer(HID_T) :: fid, gid
    integer(HSIZE_T), dimension(2) :: dims

    integer :: impc
    integer :: error

    integer, dimension(:,:), allocatable :: p_mask, f_mask

    !!Initialize HDF5 interface
    call h5open_f(error)

    !!Open HDF5 file
    call h5fopen_f(inputs%equilibrium_file, H5F_ACC_RDONLY_F, fid, error)

    !!Open PLASMA group
    call h5gopen_f(fid, "/plasma", gid, error)

    !!Read in interpolation grid
    call h5ltread_dataset_int_scalar_f(gid, "/plasma/nr", inter_grid%nr, error)
    call h5ltread_dataset_int_scalar_f(gid, "/plasma/nz", inter_grid%nz, error)

    allocate(inter_grid%r(inter_grid%nr),inter_grid%z(inter_grid%nz))
    allocate(inter_grid%r2d(inter_grid%nr,inter_grid%nz))
    allocate(inter_grid%z2d(inter_grid%nr,inter_grid%nz))
    allocate(p_mask(inter_grid%nr,inter_grid%nz))
    allocate(f_mask(inter_grid%nr,inter_grid%nz))

    dims = [inter_grid%nr, inter_grid%nz]
    call h5ltread_dataset_double_f(gid, "/plasma/r", inter_grid%r, dims(1:1), error)
    call h5ltread_dataset_double_f(gid, "/plasma/z", inter_grid%z, dims(2:2), error)
    call h5ltread_dataset_double_f(gid, "/plasma/r2d", inter_grid%r2d, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/z2d", inter_grid%z2d, dims, error)

    inter_grid%dr = abs(inter_grid%r(2)-inter_grid%r(1))
    inter_grid%dz = abs(inter_grid%z(2)-inter_grid%z(1))
    inter_grid%da = inter_grid%dr*inter_grid%dz

    if(inputs%verbose.ge.1) then
        write(*,'(a)') '---- Interpolation grid settings ----'
        write(*,'(T2,"Nr: ",i3)') inter_grid%nr
        write(*,'(T2,"Nz: ",i3)') inter_grid%nz
        write(*,'(T2,"dA: ", f5.2," [cm^2]")') inter_grid%da
        write(*,*) ''
    endif

    !!Read in plasma parameters
    allocate(equil%plasma(inter_grid%nr,inter_grid%nz))

    call h5ltread_dataset_double_f(gid, "/plasma/dene", equil%plasma%dene, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/te", equil%plasma%te, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/ti", equil%plasma%ti, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/zeff", equil%plasma%zeff, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/vr", equil%plasma%vr, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/vt", equil%plasma%vt, dims, error)
    call h5ltread_dataset_double_f(gid, "/plasma/vz", equil%plasma%vz, dims, error)
    call h5ltread_dataset_int_f(gid, "/plasma/mask", p_mask, dims,error)

    impc = inputs%impurity_charge
    where(equil%plasma%zeff.lt.1.0)
        equil%plasma%zeff = 1
    endwhere

    where(equil%plasma%zeff.gt.impc)
        equil%plasma%zeff = impc
    endwhere

    where(equil%plasma%dene.lt.0.0)
        equil%plasma%dene = 0.0
    endwhere

    where(equil%plasma%te.lt.0.0)
        equil%plasma%te = 0.0
    endwhere

    where(equil%plasma%ti.lt.0.0)
        equil%plasma%ti = 0.0
    endwhere

    equil%plasma%denimp = ((equil%plasma%zeff-1.d0)/(impc*(impc-1.d0)))*equil%plasma%dene
    equil%plasma%denp = equil%plasma%dene - impc*equil%plasma%denimp

    !!Close PLASMA group
    call h5gclose_f(gid, error)

    !!Open FIELDS group
    call h5gopen_f(fid, "/fields", gid, error)

    allocate(equil%fields(inter_grid%nr,inter_grid%nz))

    !!Read in electromagnetic fields
    call h5ltread_dataset_double_f(gid, "/fields/br", equil%fields%br, dims, error)
    call h5ltread_dataset_double_f(gid, "/fields/bt", equil%fields%bt, dims, error)
    call h5ltread_dataset_double_f(gid, "/fields/bz", equil%fields%bz, dims, error)
    call h5ltread_dataset_double_f(gid, "/fields/er", equil%fields%er, dims, error)
    call h5ltread_dataset_double_f(gid, "/fields/et", equil%fields%et, dims, error)
    call h5ltread_dataset_double_f(gid, "/fields/ez", equil%fields%ez, dims, error)
    call h5ltread_dataset_int_f(gid, "/fields/mask", f_mask, dims,error)

    !!Calculate B field derivatives
    call deriv(inter_grid%r, inter_grid%z, equil%fields%br, equil%fields%dbr_dr, equil%fields%dbr_dz)
    call deriv(inter_grid%r, inter_grid%z, equil%fields%bt, equil%fields%dbt_dr, equil%fields%dbt_dz)
    call deriv(inter_grid%r, inter_grid%z, equil%fields%bz, equil%fields%dbz_dr, equil%fields%dbz_dz)

    !!Close FIELDS group
    call h5gclose_f(gid, error)

    !!Close file
    call h5fclose_f(fid, error)

    !!Close HDF5 interface
    call h5close_f(error)

    allocate(equil%mask(inter_grid%nr,inter_grid%nz))
    equil%mask = 0.d0
    where ((p_mask.eq.1).and.(f_mask.eq.1)) equil%mask = 1.d0
    if (sum(equil%mask).le.0.d0) then
        write(*,'(a)') "READ_EQUILIBRIUM: Plasma and/or fields are not well defined anywhere"
        stop
    endif

end subroutine read_equilibrium