fieldmat_type Derived Type

type, private :: fieldmat_type

This is the top level object, consisting of a collection of ky blocks. Currently only support: 0 : Invert the field matrix In the future may wish to do things like LU decomposition etc. Could also do things like eigenvalue analysis etc.


Contents

Source Code


Components

Type Visibility Attributes Name Initial
type(ky_type), public, dimension(:), allocatable :: kyb

The ky blocks

integer, public :: naky

Number of ky blocks

integer, public :: ntheta0
integer, public :: npts

Total number of theta grid points on all extended domains

integer, public :: nbound

Number of ignored boundary points

logical, public :: is_local

Does this supercell have any data on this proc?

logical, public :: is_empty

Have we got any data for the fieldmat on this proc?

logical, public :: is_all_local

Is all of this supercells data on this proc?

type(comm_type), public :: fm_sub_all

Sub communicator involving all processors with fieldmat

type(comm_type), public :: fm_sub_headsc_p0

Sub communicator involving the supercell heads and proc0

integer, public :: prepare_type = 0

What sort of preparation do we do to the field matrix

integer, public, dimension(:,:), allocatable :: heads

An array that holds the iproc (from comm world) of the supercell head for a given ik,it point


Type-Bound Procedures

procedure, private, :: deallocate => fm_deallocate

  • private subroutine fm_deallocate(self)

    Deallocate storage space

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: allocate => fm_allocate

  • private subroutine fm_allocate(self)

    Allocate storage space

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: debug_print => fm_debug_print

  • private subroutine fm_debug_print(self)

    Debug printing

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(in) :: self

procedure, private, :: get_field_update => fm_get_field_update

  • private subroutine fm_get_field_update(self, phi, apar, bpar)

    A routine to calculate the update to the fields

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self
    complex, intent(inout), dimension(:,:,:) :: phi
    complex, intent(inout), dimension(:,:,:) :: apar
    complex, intent(inout), dimension(:,:,:) :: bpar

procedure, private, :: init => fm_init

  • private subroutine fm_init(self)

    Initialise the field objects

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: reset => fm_reset

  • private subroutine fm_reset(self)

    A routine to reset the object

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: populate => fm_populate

  • private subroutine fm_populate(self)

    Find the response of g to delta-fn field perturbations and store at the row level

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: init_next_field_points => fm_init_next_field_points

  • private subroutine fm_init_next_field_points(self, field, pts_remain, kwork_filter, ifl)

    Initialise the next set of delta functions, find the response and store in the appropriate rows.

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self
    complex, intent(inout), dimension(-ntgrid:ntgrid,ntheta0,naky) :: field
    integer, intent(inout) :: pts_remain
    logical, intent(inout), dimension(ntheta0, naky) :: kwork_filter
    integer, intent(in) :: ifl

procedure, private, :: prepare => fm_prepare

  • private subroutine fm_prepare(self)

    Prepare the field matrix for calculating field updates

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: make_subcom_1 => fm_make_subcom_1

  • private subroutine fm_make_subcom_1(self)

    Create all the necessary subcommunicators

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: make_subcom_2 => fm_make_subcom_2

  • private subroutine fm_make_subcom_2(self)

    Create the secondary subcommunicators

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: calc_sc_heads => fm_calc_sc_heads

  • private subroutine fm_calc_sc_heads(self)

    Record the process that is the supercell head associated with a given ik,it point Supercells have a single ik but may have many it points so a given supercell head may be responsible for a range of points

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: gather_fields => fm_gather_fields

  • private subroutine fm_gather_fields(self, gf_ph, gf_ap, gf_bp, ph, ap, bp)

    Redistribute fields from gf_lo to g_lo. This first involves unpacking the calculated fields from tmpsum on the supercell heads to the field arrays and then sending them back to the processes that calculated them in gf_lo and the processes that need them in g_lo. DD>TAGGED:Worth looking at improving this bad memory pattern

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_bp
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: bp

procedure, private, :: gather_init_fields => fm_gather_init_fields

  • private subroutine fm_gather_init_fields(self, gf_phnew, gf_apnew, gf_bpnew, gf_ph, gf_ap, gf_bp, phnew, apnew, bpnew, ph, ap, bp)

    Redistribute initial fields from gf_lo to g_lo. These are fields that have not been reduced to the supercell heads, just calculated on each process that owns gf_lo points

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_phnew
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_apnew
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_bpnew
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_bp
    complex, intent(inout), dimension(-ntgrid:,:,:) :: phnew
    complex, intent(inout), dimension(-ntgrid:,:,:) :: apnew
    complex, intent(inout), dimension(-ntgrid:,:,:) :: bpnew
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: bp

procedure, private, :: scatter_fields => fm_scatter_fields

  • private subroutine fm_scatter_fields(self, gf_ph, gf_ap, gf_bp, ph, ap, bp)

    Redistribute the fields from g_lo to gf_lo This is used to take the initial fields that are calculated in initialisation using the fields_implicit functionality, and therefore in g_lo, and distribute them to gf_lo layout as well so that we can start the advance_gf_local routine with the correct field data in both formats.

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_bp
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: bp

procedure, private, :: unpack_to_field => fm_unpack_to_field

  • private subroutine fm_unpack_to_field(self, ph, ap, bp)

    A routine to unpack the supercell tmp_sum vectors to full field arrays

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ph
    complex, intent(inout), dimension(-ntgrid:,:,:) :: ap
    complex, intent(inout), dimension(-ntgrid:,:,:) :: bp

procedure, private, :: write_debug_data => fm_write_debug_data

  • private subroutine fm_write_debug_data(self)

    Write some debug data to file

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(in) :: self

procedure, private, :: set_is_local => fm_set_is_local

  • private subroutine fm_set_is_local(self)

    Just work out the locality (also sets is_empty etc. but is not intended for this)

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: count_subcom => fm_count_subcom

  • private subroutine fm_count_subcom(self)

    Count how many subcommunicators will be created

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(inout) :: self

procedure, private, :: getfieldeq_nogath => fm_getfieldeq_nogath

  • private subroutine fm_getfieldeq_nogath(self, phi, apar, bpar, fieldeq, fieldeqa, fieldeqp)

    FIXME : Add documentation

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(in) :: self
    complex, intent(in), dimension (-ntgrid:,:,:) :: phi
    complex, intent(in), dimension (-ntgrid:,:,:) :: apar
    complex, intent(in), dimension (-ntgrid:,:,:) :: bpar
    complex, intent(inout), dimension (-ntgrid:,:,:) :: fieldeq
    complex, intent(inout), dimension (-ntgrid:,:,:) :: fieldeqa
    complex, intent(inout), dimension (-ntgrid:,:,:) :: fieldeqp

procedure, private, :: getfieldeq1_nogath => fm_getfieldeq1_nogath

  • private subroutine fm_getfieldeq1_nogath(self, phi, apar, bpar, antot, antota, antotp, fieldeq, fieldeqa, fieldeqp)

    FIXME : Add documentation

    Arguments

    Type IntentOptional Attributes Name
    class(fieldmat_type), intent(in) :: self
    complex, intent(in), dimension (-ntgrid:,:,:) :: phi
    complex, intent(in), dimension (-ntgrid:,:,:) :: apar
    complex, intent(in), dimension (-ntgrid:,:,:) :: bpar
    complex, intent(in), dimension (-ntgrid:,:,:) :: antot
    complex, intent(in), dimension (-ntgrid:,:,:) :: antota
    complex, intent(in), dimension (-ntgrid:,:,:) :: antotp
    complex, intent(inout), dimension (-ntgrid:,:,:) :: fieldeq
    complex, intent(inout), dimension (-ntgrid:,:,:) :: fieldeqa
    complex, intent(inout), dimension (-ntgrid:,:,:) :: fieldeqp

procedure, private, nopass :: update_fields => fm_update_fields

  • private subroutine fm_update_fields(phi, apar, bpar, gf_phi, gf_apar, gf_bpar)

    Update the fields using calculated update We update both the g_lo (phi,apar,bpar) and gf_lo (gf_phi,gf_apar,gf_bpar) fields as these are needed to make sure we can progress the fields in g_lo and gf_lo both times and thereby not have to transfer the full fields for all processes every iteration, only the required updates. With some refactoring it should be possible to do all the work on a single set of field arrays as they are the same size regardless of the layout being used.

    Arguments

    Type IntentOptional Attributes Name
    complex, intent(inout), dimension(-ntgrid:,:,:) :: phi
    complex, intent(inout), dimension(-ntgrid:,:,:) :: apar
    complex, intent(inout), dimension(-ntgrid:,:,:) :: bpar
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_phi
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_apar
    complex, intent(inout), dimension(-ntgrid:,:,:) :: gf_bpar

procedure, private, nopass :: update_fields_newstep => fm_update_fields_newstep

Source Code

  type, private :: fieldmat_type
     type(ky_type), dimension(:), allocatable :: kyb !< The ky blocks
     integer :: naky !< Number of ky blocks
     integer :: ntheta0 
     integer :: npts !< Total number of theta grid points on all extended domains
     integer :: nbound !< Number of ignored boundary points
     logical :: is_local !< Does this supercell have any data on this proc?
     logical :: is_empty !< Have we got any data for the fieldmat on this proc?
     logical :: is_all_local !< Is all of this supercells data on this proc?
     type(comm_type) :: fm_sub_all !< Sub communicator involving all processors with fieldmat
     type(comm_type) :: fm_sub_headsc_p0 !< Sub communicator involving the supercell heads and proc0
     integer :: prepare_type=0 !< What sort of preparation do we do to the field matrix
     integer, dimension(:,:), allocatable :: heads !<  An array that holds the iproc (from comm world) of the supercell head for a given ik,it point
     !! Currently only support:
     !! 0   : Invert the field matrix
     !! In the future may wish to do things like LU decomposition etc.
     !! Could also do things like eigenvalue analysis etc.
   contains
     private
     procedure :: deallocate => fm_deallocate
     procedure :: allocate => fm_allocate
     procedure :: debug_print => fm_debug_print
     procedure :: get_field_update => fm_get_field_update
     procedure :: init => fm_init
     procedure :: reset => fm_reset
     procedure :: populate => fm_populate
     procedure :: init_next_field_points => fm_init_next_field_points
     procedure :: prepare => fm_prepare
     procedure :: make_subcom_1 => fm_make_subcom_1
     procedure :: make_subcom_2 => fm_make_subcom_2
     procedure :: calc_sc_heads => fm_calc_sc_heads
     procedure :: gather_fields => fm_gather_fields
     procedure :: gather_init_fields => fm_gather_init_fields
     procedure :: scatter_fields => fm_scatter_fields
     procedure :: unpack_to_field => fm_unpack_to_field
     procedure :: write_debug_data => fm_write_debug_data
     procedure :: set_is_local => fm_set_is_local
     procedure :: count_subcom => fm_count_subcom
     procedure :: getfieldeq_nogath => fm_getfieldeq_nogath
     procedure :: getfieldeq1_nogath => fm_getfieldeq1_nogath
     procedure, nopass :: update_fields => fm_update_fields
     procedure, nopass :: update_fields_newstep => fm_update_fields_newstep
  end type fieldmat_type